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Chapter 1 

Introduction 



Most stars lose mass. But the importance of the stellar mass loss for the fate of the star and the rest of the 
universe differs quite from star to star. Our sun have a negligible mass los rate while the mass loss from O, B, 
and Wolf-Rayet stars is very important. 

O and B stars are known for a long time. While the first Wolf-Rayet stars were discovered by Wolf & Rayet 
not before 1867 [7C(] . The reason for that is that they are so few. Only about 160 Wolf-Rayet stars are known in 
our galaxy today. Some more we know in the Magellanic clouds and in M33. The most prominent spectroscopical 
feature of Wolf-Rayet stars are their strong emission lines, while the spectra of most normal stars like our sun are 
dominated by absorption lines. Similar spectra can be seen in other, less luminous objects like planetary nebula. 
The common reason for the strong emission lines is a strong, optically thick outflow of matter. In this thesis 
we consider only classical Wolf-Rayet stars from the population I when we use the terminus "Wolf-Rayet star" . 
Although planetary nebulae have a strong mass outflow as well, the underlying physics and the stellar evolution 
are different. Wolf-Rayet stars play an especially interesting role in astrophysics. They have the strongest stellar 
winds. And their winds are strongly enhanced with heavy elements. Their winds contain nearly no hydrogen 
anymore. The Wolf-Rayet stars of the WN subclass show the chemical composition of the products of the CNO 
cycle in their winds. While the stars of the WC subclass already show the products of helium burning in their 
winds. In both cases we have the bare nuclei of more massive progenitor stars. But massive stars do not only 
enrich the interstellar medium with heavy elements. Additionally the stellar winds carry momentum and kinetic 
energy. These are very important to trigger e.g. star formation. In this thesis we will only consider blue stars, 
like O, B, and Wolf-Rayet stars when we talk about "massive stars." There are massive stars in the red part of 
the Hertzsprung-Russel diagram as well. But their winds contain dust and are driven in a completely different 
way. Therefore they are excluded in this thesis. 

Significant mass loss as a common phenomenon in massive stars is known only for a short time, because only 
few massive stars like Wolf-Rayet stars show the typical P-Cygni line profile which is the spectral indicator for 
strong mass loss. A P-Cygni line is the combination of an emission line with an absorption line connected to the 
blue wing of the emission line. Both lines come from the same atomic transition. The absorption is caused in 
the line of sight by the wind material moving towards the observer. This causes the blue-shift of the absorption 
line. From the blue edge of the absorption line we can therefore derive the maximum velocity in the outflowing 
matter. The emission line is produced in the wind material which leaves the star roughly perpendicular to the 
line of sight. Only few, weak P-Cygni lines lie in the optical range. Therefore only few stars with very strong 
winds like Wolf-Rayet stars are known for long time to have significant mass loss. Thus significant mass loss 
was considered to be a rare phenomenon in stars. In 1967 Morton used rockets to observe massive stars 

in the ultraviolet (UV) range. He found that many massive stars show P-Cygni lines in the UV-range. And so 
mass loss became known as a phenomenon common in massive stars. It became also clear that the absorption of 
stellar light in atomic lines plays the dominant role in driving these winds. Lucy & Solomon |39| did pioneering 
theoretical work in this field. Based on their results Castor, Abbott & Klein jl6| developed a genius model to 
describe the line driving mechanism in a simple but in most cases accurate way. We will make extensive use of 
their model in this thesis. 

Today we have much more detailed observations of massive stars in many frequency ranges. This gives us 
a modern and more detailed picture of their winds. The wind material can today be detected in the radio 
range, as well. This is more than just a confirmation of the results from P-Cygni lines because normal wind 
material produces only thermal radio emission. But many massive stars show nonthermal radio emission ^, [|. 
This is synchrotron emission, which is a clear evidence for magnetic fields in these winds. Additionally an 
acceleration mechanism is required to produce high energy particles for the synchrotron emission. The most 
obvious mechanism is Fermi acceleration in shocks. Shocks appear where the wind collides with the interstellar 
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medium or with the wind of a neighbouring star (e.g. in binaries). Additionally shocks can be produced by the 
instability of the wind driving by radiation. Further evidence for the existence of shocks in these winds come 
from the observation of X-rays in these winds [^6| . Detailed analysis of the stellar spectra show strong evidence 
for wind variability as well. Thus we can not describe the wind of massive stars as stationary, smooth outflows 
any more. Modern wind models should describe in more detail how these perturbations of the smooth wind are 
created, how they are amplified, and which consequences they have for the overall wind. From these questions 
we can learn a lot about these winds. 

The first observations of stellar winds were made for our own sun, where due to the short distance we can 
study wind phenomena much easier and in far more detail than for other stars. L. Bicrmann was the first who 
recognized that the tails of the comets, which always point away from the sun, are a clear indication for a 
permanent matter outflow from the sun |^|, Q . Based on this information Parker developed the first quantitative 
model for the solar wind. He showed that the hot corona of the sun (~ 10 6 K) must have a permanent matter 
outflow, because gravity and the interstellar medium can not compensate the high thermal pressure. During the 
years this basic model has been improved by theory and observation in many ways. Today solar physics is a 
broad field in astrophysics. Most of these improvements are too detailed to be used for winds from other stars. 
But one important model, first developed for our sun, will be used in this thesis extensively. It is the model of 
Weber & Davis |39j for the solar magnetic field. In the solar wind the magnetic field plays only a subordinate 
role. It merely controls the solar spin-down. The key message of this thesis is that the magnetic field plays an 
much more important role in the wind physics of massive stars. Massive stars rotate faster and we have many 
indirect evidence that they have strong magnetic fields ||. This can enhance the wind drastically. While it 
seems possible to describe the winds of O and B stars with a purely radiative model [[32, 59 , this attempt failed 
for Wolf-Rayet stars. 

Models which use the magnetic field as additional ingredient for the winds of massive stars have been devel- 
oped earlier |2C], [jl, 62 1. But these models suffered from various problems. The strongest limitation to magnetic 
wind models is the so called spin-down problem. The basic idea of the Weber & Davis models is that open field 
lines are fixed in the stellar envelope and extend outwards. These field lines co-rotate with the star. Through 
the Lorentz force they help to push the wind material outside. Through this process the star will lose an extra 
amount of angular momentum. If the magnetic field is to strong the star will spin down very rapidly. In this 
case the Weber & Davis mechanism can not play a role in further stages of the stars evolution. Therefore is 
is important to keep the evolution of a star though its whole life in mind. It is e.g. not clear whether massive 
stars have a magnetic field through their whole lifetime. It is possible that young, massive stars on the main 
sequence have a purely radiation driven wind, and therefore radiative wind models are so successful for these 
stars. Evolved stars like Wolf-Rayet stars, whose wind are not explained yet, may have a magnetic field, which 
was either generated in a later stage of the stellar evolution or was hidden under the outer layers of the stellar 
envelope blown away now. These are very speculative ideas. But we are still in an early stage of understanding 
massive stars. 

Another important aspect of magnetic fields in winds of massive stars is the production of cosmic rays. We 
know from the observation of nonthermal radio emission that we have shocks and a magnetic field in the winds. 
If the star rotates the magnetic field will be bent into azimuthal direction far away from the star. If we combine 
such a field configuration with shock running in radial direction we get an optimal configuration for particle 
acceleration due to the Fermi mechanism. The high energy particles produced here gyrate in the magnetic field 
and produce the observed nonthermal radio emission. But high energy particles, especially protons and nucleons, 
can also escape from the system and contribute to the galactic cosmic ray spectrum. This is the basic concept 
of Biermann & Cassinelli g, ^| , who combined several aspects of magnetic fields in the winds of massive stars 
into a self-consistent picture. This work gave rise to the initial concept of this thesis. 

In this thesis we analyze several aspects of magnetic fields in stellar winds. The fundament for these models 
are the improved version of the fast magnetic rotator model developed by Biermann & Cassinelli || ^) . The first 
major part of this thesis is a model for the influence of the magnetic field on the radiative instability in the wind. 
Already very early it became clear that line driven winds are unstable |39], ^6) . But all previous attempts to use 
these instabilities to explain the thick and fast winds, we observe in massive stars, failed. We found that this is 
due to the neglect of the magnetic field. Previous work concentrated on the effect of radiation and found that 
the instability produces inward running waves p2[ . These waves can help to explain many details in the spectra 
of massive stars. But they do not help to explain basic wind parameters like the high observed terminal wind 
velocities. We found that a magnetic field turns the waves into outward running waves |6^] , which can explain 
the high observed terminal velocities with moderate magnetic fields. The latter is very important to keep the 
spin-down problem under control. 

In the second major part of this thesis we extend the basic wind model to take the wind outside the equatorial 
plane into account. Older wind models separate into two classes. The first class of models ignores rotation and 
the magnetic field. These models treat all quantities as only radius dependent. Thus they are one dimensional 
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and strictly spherical. Most purely radiation driven wind models belong to this class. The second class of models 
treats rotation and often magnetic fields. They assume only rotational symmetry but restrict themselves to the 
equatorial plane. Through this restriction they are essentially one dimensional as well. Strong assumptions are 
necessary to estimate the wind outside the equatorial plane. Our wind model keeps the simplicity of an one 
dimensional wind model but allows to explore the wind outside the equatorial plane. This gives us much more 
confidence in wind models with rotation. The second important aspect of our extension is, that this model can 
now be used to calculate winds which do not blow in the equatorial plane at all. For such winds the old equatorial 
models are meaningless. This is the large family of winds from accretion disks, which often turn into jets close 
to the polar axis. Such winds appear in many astrophysical objects from young stellar objects (YSO) to active 
galactic nuclei (AGN). Due to this broad applicability we expect many interesting results from our wind model 
in the future. 

In Chap. |^ we start our discussion with some basic concepts explained at the wind model of Parker. In 
Chap. |^ we give a short introduction to the theory of line driven winds. Since this thesis concentrates on the role 
of the magnetic field only the basic concepts for our own usage are given here. In Chap. ^ the wind equations 
for the equatorial plane are derived. These equations are already slightly generalized compared to the equations 
of Biermann & Cassinelli by the fact that they incorporate the effect of a wind compressed or diluted in the 
equatorial plane. In Chap. [3| we combine the results from the two previous chapters for some initial wind models, 
which allow the comparison with previous models. In Chap. ^ we derive our model for waves in magnetic winds 
from massive stars and give numerical results for a generic model star. Chap. (?] contains our model for the 
wind outside the equatorial plane. Some initial numerical results are given as well. In Chap. ^| we summarize 
our results in this thesis. The appendices contain some extra numerical and mathematical material for the 
non-equatorial wind model from Chap. 



Chapter 2 

Introduction to stellar wind theory 



The first stellar wind models were developed for our sun. L. Biermann j^, [?[] first proposed a steady hydro- 
dynamical outflow from the sun to explain the direction of the cometary tails. Based on this idea Parker |5S| ] 
showed that the sun indeed can not have a stationary atmosphere. Due to the high temperature of the corona 
a stationary atmosphere would have a high finite pressure at infinity. This pressure could not be compensated 
by the interstellar medium. A continuous expansion of the solar atmosphere is the direct consequence. Parker's 
solar wind is only driven by thermal pressure. This leads to a simple model which nevertheless gives important 
insight in the physics of stellar winds in general. Therefore the important parts of Parker's model will be outlined 
in this chapter. 

The solar wind has a velocity of approximately 300 — 800km/s at the earth's orbit. Even very early type 
stars produce only a wind of up to 3500 km/s much smaller than the speed of light. Therefore nonrelativistic 
physics is sufficient to describe these stellar winds and will be used in this thesis. 

For this introducing discussion to stellar winds it is sufficient to study the average, large scale structure of 
the solar wind. This allows us to use the theory of hydrodynamics, and later magnetohydrodynamics, as basis 
for our models. Hydrodynamics can be used whenever the relevant scales of the model are larger than the mean 
free path A for collisions between molecules. At the base of the solar wind we have roughly A(3i? Q ) = 4- 10~ 2 i?©. 
At the earth's orbit A has increased to A(1AU) = 1.5AU. Here and further out rigorous treatment using the 
Fokker-Planck equation would be better. But this is beyond the scope of this thesis. This problem is relaxed by 
the fact that the solar wind has a magnetic field. The magnetic field forces the ionized wind particles to gyrate 
around the field lines. This effect increases the probability to hit an other particle and allows the particles to be 
scattered at field irregularities, so that the hydrodynamical equilibrium is reestablished. Furthermore the wind 
of early type stars is much denser than the solar wind. This leads much shorter mean free paths and justifies 
the usage of hydrodynamics as well. 

The hydrodynamical wind is described by the conservation of mass 



dp 
dt 



V-(pv) = 0, (2.1) 



momentum 



and energy 



/9 |^ + p(vV)v = F ext , (2.2) 

p (ji + (v ' v)e ) + P(v ' v) = ~ Vqc ' (2 - 3) 

where p is the matter density, t the time, v the velocity of the bulk matter, F oxt the sum of all external force 
densities, e the internal energy per unit mass, P the thermal pressure, q c = -k c VT the conductive energy flux, 
k c the thermal conductivity, and T the temperature. If we assume that the star does not rotate and that the 
wind is stationary, the whole system reduces to one dimension 

M = 4npv r r 2 = const (2.4) 

pv r v' r = Fext, (2.5) 

where ' denotes the derivative with respect to the radius r. M is the mass loss rate and v r is the radial 
component of the wind velocity v. The external forces are created by radiation pressure from the stellar light, 
thermal pressure, interaction of the gas with waves and shocks, and gravitation. We will discuss radiation 
pressure as a major component of early type star wind in Chap. [| The creation of waves and shocks by the 
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instability of radiation pressure in a magnetic wind is the subject of Chap. [| Other sources for sho cks and waves 



are e.g. discussed in |lq| . The influence of waves and shocks on the average wind through Eq. 2.2 was analysed 
by Koninx ]3C| | . The gravitational potential derives from Poissons equation 

V 2 $ = -AirGp, (2.6) 

where G is the gravitational constant. The contribution of the wind mass 

pdV^M — (2.7) 

v r 

is negligible compared to the stellar mass, where is the radius of the heliosphere. For the sun we find 
-Wwind ss 1O _12 M0, where M is the solar mass. All normal stars are also reasonably spherically symmetric. We 
can therefore approximate the gravitational potential by 

• (2.8) 
r 

where M is the stellar mass. We now express the thermal pressure P by the isothermal speed of sound v s 

P = v 2 sP (2.9) 
< = ^ (2-10) 



where m is the average mass of the wind particles. If we use gravitation and thermal pressure as the only relevant 
forces, as Parker did 

Fcxt = -p*' - P', (2.11) 

we can write the momentum equation as 

f(r ) 2£ _ y» - QM 
Vr v' r = -l^-T = ~ S 2 ^ ■ (2.12) 

g{r,v r ) 1-% 

v r 

In the general case this ordinary differential equation has to be solved simultaneously with an equation for T(r) 
which can be derived from Eq. |2.3|. This can be found in pq, |45|. For our discussion it is sufficient to analyze 



the isothermal approximation. In this case g in Eq. 2.12 depends only on v r . For v r > v s , g is larger than zero 

(2.13) 



and otherwise smaller than zero. f(r) changes its sign at the Parker radius 

GM 



This leads to a classical family of solutions. The members with v r > are sketched in Fig. 2.1, If a solution 
reaches \v r \ = v s at a radius r ^ r c , v' r becomes infinite. This leads to double valued solutions like solutions 
5&6. These are unphysical for a stationary wind. If a solution passes through the critical point |u r | = v s , r = r c , 
l'Hospital's rule gives 

«;m= ± ||, (2.i4) 

which leads to two solutions. Solution 1 is subsonic at the solar surface and continously accelerates up to 
supersonic velocities. This solution was proposed by Parker for the solar wind. Solution 2 describes a wind 
which has for some reason a very high velocity at the solar surface. This does not fit the observation there. The 
same is true for solution 4, which is supersonic at all radii. Solution 3 did indeed fit all observations in 1958. 
It was proposed by Chamberlain |l7j as solar breeze. Finally the question was solved in favour of the transonic 
winds by observations from satellites, which found a supersonic plasma at the earth's orbit. The solutions with 
negative velocities can be interpreted as accretion scenarios [^9[ Q . 

Critical points are very useful for the solution of wind equation, because they contain extra information. 



Equation 2.12 is a single ordinary differential equation of first order. For such an equation it is normally 
necessary to specify v r at one point in order to select a unique solution from the infinite family of possible 
solutions. Once wc know from physical arguments that the solution has to go through the critical point, only 
a small finite number of solutions remain. This becomes important in the case of stellar winds, where less 
parameters are accessible for direct observation than in the case of our sun. 

The isothermal model has the advantage that neglecting the energy balance simplifies the mathematics. 
Physically it means that we implicitly introduce additional energy sources and sinks to maintain the constant 
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Figure 2.1: The solutions of Parker's isothermal wind equation 



temperature. Close to the sun, where the temperature changes only by one order of magnitude between the 
corona and the earth's orbit, this is a reasonable approximation. When we go further out the wind has to cool 



down adiabatically finally. In the isothermal case Eq. 2.12 can be readily integrated to 



2 

w r \ _ . / v r \ . . / r \ r. 



T.) - 21 "UJ= 41 "UJ +4 - + a (2 - 15) 

For r — > oo and v r 3> v s this leads to v r = 2v s In ' 5 (r/r c ) — > oo, because we add an infinite amount of thermal 
energy to maintain the constant temperature at infinity. 

Real stellar winds do not extend to infinity, but they collide with the interstellar medium. This happens 
when the ram pressure equals the pressure of the interstellar medium. Therefore all numerical models for stellar 
wind compute the wind only up to a large but finite radius. This radius should be chosen so that it is larger 
than any radius relevant in the wind. In Chap. ^ we will use u ~ r _1 as spatial coordinate. In this case it is 
mathematically possible to extend the computation to infinity in the mathematical sense. This can not be done 
using the isothermal approximation. But since the wind is terminated in the real world by the collision with 
the interstellar medium anyway, we will compute our wind models only to a large radius using the isothermal 
approximation. 

In this chapter we established the fundamental framework for our stellar wind models. We will use nonrela- 
tivistic hydrodynamics or magnetohydrodynamics (MHD) in the stationary or quasi stationary form to describe 
stellar winds. This thesis emphasizes the role of the magnetic field. Thus we treat the wind as isothermal and 
do not go into the details of radiation transfer theory. 



Chapter 3 

Stellar radiation and winds 



3.1 Introduction 

Photons carry momentum and energy. Both can be exchanged with matter. Therefore it seems quite obvious 
that the stellar light can play an important role in stellar wind physics. Unfortunately, the detailed treatment 
of this process is very complicated. 

We start our discussion with some basic definitions from radiation theory following the description of Mihalas 
PHl and Falcke |jTs| . The specific radiation intensity 

v dvdtn- da<i£! 

is the differential energy dE ra & per unit frequency v and time t, which is emitted from an area of size da oriented 
in direction da into an solid angle dCl in direction n. If the radiation with solid angle angle dfl and direction n 
passes through the area da along a distance ds, and the volume dV = dads is filled with matter of density p, 
the energy 

dEabs = I v pK v d£l dv dt n • da ds (3-2) 
is absorbed. The opacity k v can be equivalently expressed by the total cross section for absorption 

a v = mK u . (3-3) 
The momentum transfer dp — n dp per unit volume and due to the absorption can be expressed by 

This assumes that all photons which are emitted from the matter, including scattered photons, are emitted 
isotropically. Therefore their net momentum is zero and can be neglected in the momentum transfer between 
radiation and matter. This absorption approximation is sufficient for our models. More elaborate radiation 



transfer models, which include non-isotropic emission, have been calculated by Gayley et al. |23 and others. 
To obtain the total radiative acceleration we have to integrate dp over frequency and the angle of incoming 
radiation: 



1 

Trad 

p 



dp (3.5) 

l v K, v wd£ldv (3-6) 

poo 

c" 1 / F v n v dv (3.7) 



o 

a and k are complicated functions of temperature, density, chemical composition, frequency, and the intensity 
I v . Since the radiation and matter influences each other while the radiation passes through many layers of 
material, we end up with an enormous nonlocal set of equations. Solving this problem rigorously leads to good 
quantitative results for the radiation transfer but also needs an enormous effort |3^, p9[ . This is beyond the 
scope of this thesis. Therefore it is necessary to simplify the problem to a level which is reasonably close to 
physical reality but simple enough. For the treatment of hot star winds at this level there are two important 
models for the interaction of radiation and wind: Thompson scattering and line absorption. 



In Sect. 3.2 we discuss Thompson scattering as a very simple mechanism for the interaction between matter 



and radiation. In Sect. 3.3 we discuss how photon absorption in atomic lines can be described in a simple model. 



And finally we discuss in Sect. 3.4 the critical point in solutions for these radiation driven winds. 



8 CHAPTER 3. STELLAR RADIATION AND WINDS 

3.2 The Thompson wind 

Due to the high temperature in the wind of hot stars and due to the strong ultraviolet radiation field, the matter 
in the wind is highly ionized. This leads to a large number of free electrons. Free electrons are, compared to 
atoms and molecules, very simple objects. Their total cross section for scattering of radiation is 

<XTh = Y r oe = 6 - 65 10~ 25 cm 2 (3.8) 
e 2 

r cc = J ( 3 - 9 ) 



where r cc is the classical electron radius. The angular dependence of the differential cross section and the 
backscattering of the photons can be neglected for our treatment (absorption approximation). To calculate the 
opacity we need to know the number of free electrons per unit mass. Since the wind can be assumed to be fully 
ionized this depends only on the ratio of protons (and electrons) to nucleons. For hydrogen this is 1, for helium 
and heavier elements it is 0.5. This leads to 

K Th = (1 + X)£2L« 0.2(1 + X) — , (3.10) 
2m p g 



where X is the number fraction of hydrogen in the wind. Since KTh is independent of z/, we can integrate Eq. 3.7 



directly. If we assume that the wind is optically thin, so that the flux F„ depends only on the geometry, we get 

■9Th = "j 2~> ( 3 - U ) 

47rr z c 

where 

L= I / F„dVldv = 4irr 2 F (3.12) 



is the total luminosity of the star. Here we assumed that the star is point like. This simplifies the integration 
over the angle of the incoming photons. Models which take the finite size of the star into account had been 
calculated for the C AK-model p9l p0| . The Thompson radiation pressure has the same radial dependence than 



gravity. So if we add Thompson scattering to the external force balance (Eq. 2.11 ) we get 



-Fext = p{9grsv + 9Tb) + -^pressure (3.13) 

GMp -{l-T)-P' (3.14) 



r _ L E ~ AttGMc (3 ' 15) 

r is the well known Eddington factor, which describes the modification to gravity by Thompson scattering. If T 
would exceed 1 the star would disintegrate, because the radiation pressure overcomes gravity. This is known as 
the Eddington limit. Since the Thompson model is so simple it can be integrated without extra effort in every 
wind model. 



3.3 The theory of Castor, Abbott & Klein 

In the last section we saw that a sufficient high opacity can drive any wind. But the Thompson opacity is not 
high enough to explain the strong winds we see in many hot stars. Already in the twenties Johnson |28j and 
Milne (4^, ^| proposed absorption in spectral lines as an additional mechanism to drive outflows from stars. The 
line opacity, especially in resonance lines, is much higher than the continuum opacity. But this is much more 
difficult to treat due to the complicated interaction between the radiation intensity, as function of frequency 
and radius, and the opacity depending on individual excitation levels of the atoms. If the wind is very thin 
(i.e. the radiation is not strongly attenuated - even in lines) we could used an opacity averaged in frequency 
like the Rosseland opacity to calculate a 'gray' wind model. This would just replace the low Thompson opacity 
by the higher Rosseland opacity. But in the wind of hot stars many lines become optically thick. The strong 
dependency of the line opacity and the radiation intensity on frequency requires to take the Doppler effect for 
accelerating winds into account. Sobolev |38| did pioneering work to find a simple formalism to describe this 
mechanism. Lucy and Solomon (3^] were the first who applied this idea to hot stars shortly after Morton [ 
had observed P-Cygni profiles in the UV-spectra of such stars. But their model used only a few resonance lines 
and therefore gave much smaller mass loss rates than observed. Castor Abbott, and Klein [[l3], [l(| (CAK) solved 
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3.3. THE THEORY OF CASTOR, ABBOTT & KLEIN 



9 



this problem by including a large number of lines. Nevertheless they managed to keep the theory simple by 
encapsulating the complicated atomic physics in a simple fit formula. This allows us to use this formula in this 
project. 

We do not need the full details of this theory. Therefore we follow here the derivation of Owocki ||5l| . The 
CAK-theory is a hydrodynamical theory. It assumes that the momentum which is absorbed by some heavy 
elements is distributed among the whole wind material. This distinguishes it from the early work by Johnson 
and Milne, who assumed the heavy elements escape from the star while the hydrogen remains bounded. The 
second assumption is that the temperature and any nonthermal population of the atomic levels are known in 
advance. This fixes the opacity of the wind material. The CAK-theory takes into account, that optically thick 
lines can shadow themselves from the stellar light, and that this effect can be suppressed by the Doppler-shift 
in the expanding wind. 

We start the derivation with a single line of opacity 

k v = n L (j)(x) (3.16) 
x = (3.17) 



fth = W , (3.18) 

where x is the deviation from the line frequency v\, in units of the Doppler shift due to the thermal motion Vth- 
Here kl is the line opacity, <j){x) is the normalized line profile function, c is the speed of light, /cb is the Boltzmann 
constant, and m p is the mass of the proton. The deviation of the mean mass m of the wind particles from the 
proton mass is later incorporated in the CAK parameter a ca k and fc C ak- If we neglect now the momentum of the 
emitted photons, we can express the radiative acceleration of this line resulting from a flux F v by the momentum 
of the absorbed photons 



gh(r, ^ l ,kl) ~ 5abs(r-, k l ) = / K^Mx- EiL " LVth e r| >< r ) dx, 

Jo \ vth J c c 



(3.19) 



where we have taken the attenuation of the flux F v by the optical depth 



r(x,r) = I K L p(f)0( x- ) df (3.20) 



R \ v tl 



kl p Mx)-^y- dx (3-21) 

.Mr) v' r 

into account, where R is the radius of the stellar surface (v r (R) « 0) and x — x — (v r (f)/vth) the frequency 
in the comoving frame. The line profile function <f>{x) is sharply peaked at x = and therefore only the small 
spatial region of the wind with \x — v r (f)/vth\ ^$ 1 will contribute to the integral in Eq. [3.21 , The thickness of 
this region is given by the Sobolev length 

Ls = ^ (3.22) 



with the velocity gradient of the wind 



dv r 

(3.23) 

dr 



The important contribution of Sobolev at this point was the observation that p[f) /v' r (f) is approximately constant 
over this region of the wind. Since we have a single line wind, where the same line causes the flux attenuation and 
the wind acceleration, both effects take place at approximately the same radius r. We can therefore approximate 
p(f)/v' r (r) by p(r) / v' r (r) . This allows us to draw this term out of the integral 



r(r,x) « r L / 4>{x)dx, (3.24) 

J x _ ? r (r) 
v th 

where we have introduced the Sobolev optical depth of a single line 



(3.25) 
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For convenience we shift now to the comoving frame at radius r by x = x — v r (r)/vth to solve Eq. 3.19 analytically 

-pLKL^th^L f°° 

9h(r, v l ,k l ) w 5 / r L 0(a;)x 

c T L J-oo 



exp I — tl y do; ) dai (3.26) 



(3.27) 



In the case of an optically thin line tl < 1 we find 



#thin — 5 . (6.2s) 



In the optically thick case tl > 1 we find 



<7thin F L V L , 

Sthick = = — — v r . (3.29) 

tl crp 

This expression is independent of kl because the thick line will absorb all the photons in its frequency range 
anyway. But ffthick depends on the velocity gradient v' r , because v' r describes how efficiently the line can escape 
its own shadow and therefore can absorb new photons at a different Doppler shifted frequency. 

The major improvement of Castor, Abbott, and Klein was to take many lines into account by summing 
Eq. 3.1£ for many lines. This includes the assumption that different lines do not overlap in frequency even 
taking the Doppler shift due to the wind velocity into account. For a real line list this computation can only be 
done numerically. It is then possible to encapsulate the atomic physics of the line list in a function M(t) called 
the force multiplier. CAK found that this force multiplier can then be approximated by 

M(t)=k cak t- a ^, (3.30) 

where. 

t = KThpVth (3.31) 

v 

is the line independent Sobolev optical depth. CAK |l(| found fc C ak = 1/30 and a ca k = 0.56. Later Abbot [g] 
tabulated improved values for fc ca k and a ca k. 

The results given in Eq. 3.30 can be obtained analytically if we assume that the lines have a flux- weighted 
number distribution that has a power law in opacity given by 

tf(*L,*L) = — =--( — ) , (3.32) 

where F is the total radiation flux and ko is a normalisation constant. We can now incorporate this expression 
in Eq. 3.19 and integrate over opacity. 

/>oo 

g ca k(r,v r ,v' r ) = / N(v L ,K L )g L (r, v L ,K L )dK L (3.33) 



(—) 

\K pVthJ 



a oak 

X 



oo poo 



oo JO 

FvthKo 



T L 



a-k-l e -T(4,r) dT h <f>(x) dx (3.34) 



(—)' 

\K pVth ) 



r(a ca k) x 



4>(x)dx) (j)(x)dx (3.35) 



Fvt^Ko r(a cak ) / v' r 



c 2 1 - a cak \K pvth 

F ZtTh 



(3.36) 



Mti), (3.37) 
c 

where T(a ca k) is the complete Gamma function. Here we find the relation between k,q and k ca k: 

fc cak = ^ ( -_-V^ k (3 . 38) 
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If we include g ca ^ and gTh in the wind equation we get different analytical conditions for a complete physical 
solution than in the case of a Parker wind. The reason for this is that g ca k depends additionally on the local 
wind acceleration. To analyze the solutions we rewrite the wind equation (Eq. 2.12) to 



A{r,v r )v' r + B{r) = C{r,v r )v' r a 



(3.39) 



D 



C 



v r 



(1 



r)-2^ 

r 



GM 

Airr 2 c \KT h pvth 
GMTk cak ( iirr 2 v, 



cak 



l 



\KThMvth) 



Equation 3.39| is an implicit differential equation, which in general can not be solved analytically, 
one of the minor numerical problems in CAK wind models. In the case of an isothermal wind we 
these equations by substituting u = — 1/r, w = v 2 ./2, and w = (dw/du): 

A(w)w + B{u) = Cw a ™ k 



(3.40) 
(3.41) 
(3.42) 

(3.43) 

But this is 
can simplify 

(3.44) 



A(w) 
B(u) 

C 



1 - ^r- 

2w 

GM{\ - r) 

GMrfe cak 



2v 2 s u 
4tt 



th 



(3.45) 
(3.46) 
(3.47) 



In the case of a cold wind without thermal pressure (v s = 0) Eq. 3.44 becomes independent of r. Therefore 
w = (dw/du) is independent of r as well. This leads to the often used velocity law 



v r (r) 



1 



I2w 
~R ' 



R 



(3.48) 
(3.49) 



where R is the stellar radius. Since in this case A = 1 and B = GM(1 — T) are larger than zero, we can have 
two (Fig. [OJ.A), one (Fig. 3.1.B), or no (Fig. 3.1.G) solutions for Co depending on the value of C(M). The case 
of one solution is special, because it has the minimum value for C(M) and therefore the maximum value of M. 
The unique value for w in this case leads to 



Voo(M n 



4nGM 

KTh^th 



acak(fccakr) 



1/acak 



_Ocak\ 



(l-a cak )/Qcak 



l-r / 



a cak 2GM{1 - Y) 



1 - a ca k 



R 



(3.50) 
(3.51) 



This case is often referred to as the analytical CAK-solution. Nevertheless this situation is not satisfying, because 
solutions with lower mass loss rates are not excluded by any physical argument. 



3.4 The CAK point 



For a unique solution wit h a unique value for M we have to take the thermal pressure into account. Since a C ak is 

c an h ave no, one, or two solutions for v' r depending on the values of A, B. The possible 
If C is zero we have the Parker or Tho mpso n wind as limiting case, 
we have a negative value for B. This requires A > as plotted in Fig. 3.1. 



between and 1, Eq. 3.1 



cases are sketched in Fig. 3.1 



At r = oo 

.D and therefore a supersonic wind. 
Close to t he s tar the gr avit y term in Eq. |3.41| will dominate and cause a positive value for B. This transition 
from Fig. 3.1.D to Fig. 3.1 A takes place in the supersonic part of the wind (^4 > 0). Otherwise we would end 
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Figure 3.1: The possible solutions of Eq. 3.39 : The left hand side of Eq. 3.39 is linear in v' r : Av' r + B, while the 
right hand side: C v'" cak d e pend s on v' r roughly like v/u£ because a ca k is between an d 1. But A and B can 
have either sign (Eqs. 3.40 & 3.41 ) depending on the values of r and v r (r). From Eq. 3.43| w e see that C is always 
positive. In Figs. A-H we sketch the left hand side and the right hand side of Eq . 3.3S (ordinate) as function 
of v' r (abscissa), r and v r are kept constant in each figure. A solution for Eq. 3.39| is giv en at every value of v' r 



where the two curves in terse ct 
And in Figs. G&H Eq. [3~39 



In Fig. A Eq. 3.3S has two solutions. In Figs. B-F Eq. 3.3S has one solution, 
has no solution. In general the coefficients A and B are not constant in a single 
wind solution. Thus for different radii different of the figures sketched here may apply in the same wind solution. 
For a physical wind solution r and v r and thus A and B are continuous everywhere. Therefore the transition 
between these different figures has to be continuous as well. A nonmagnetic wind passes from the base of the 
wind to infinity through the following sequence: C-E-A-B(here is the CAK critical point)-A-D. See Sect. 3.4 
for details. The case of Fig. F appears only if magnetic fields are included (Chap. ||). In this case A and B are 
infinite at the Alfvenic point. 
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up with the situation of Fig. [3J.H where no solution for v' r can be found. In the case of Fig. 3T A we have two 
possible solutions. But the continuous transition from Fig. |3.l| .D allows only the larger value of v' r . 

If we use the CAK-model in the subsonic part of the wind we get a negative value for A. Therefore B must 
be positive and we have the unique solution for v' r , which is sketched in Fig. 3.I.C. At the sonic point (A = 0, 
Fig. " 



3.1 



.E) the solution switches to Fig. 3.1. A. But now the proper solution for v' r is the lower value, because 
only this matches continously the unique solution at the sonic point. 

If we want a continuous steady state solution, which extends from the base of the wind (v r <C v s ) to r — oo, 
we have to connect the solution branches discussed above in a proper way. This can only happen if the solution 

Only here we can shift v' r continously from the lower value in 



passes through a point where Fig. [3.1 
Fig. OA to Fig. EO 



.B applies. 

.B and then to the higher value in Fig. 3.1. A. This point is the critical CAK-point. It fixes 
the whole solution. The best way to solve Eq. 3.39 numerically is to find the critical CAK-point and to integrate 
Eq. 3.39 from there to the base of the wind and to r = oo. 

, and v'. Therefore we need four conditions 



The critical CAK-point is specified by four quantities: M, r c 



to find the critical CAK-point. The first one is Eq. 3.39. From Fig. 3.1.B we see that both sides of Eq. [3.39 have 



the same inclination when they are plotted as function of v' r . This happens only in Fig. 
the second condition: 

A(r c ,v rc ) = a cak C(r c ,v rc )v' r ° c 



3.1.B and is therefore 



(3.52) 



A physical solution for v r should be smooth everywhere. This allows us to differentiate Eq. |3.39| totally with 
respect to r: 



OA OA , 

dr dv r 1 



dC_ 

dr 



dC 
av r 



Av'l 



,,' Q cak 



dB 

dr 



Cv'S c 



At the critical CAK-point we can use Eq. 3.52 to eliminate v 



8A 
dr 



OA 

dv r 



dB 
dr 



dB 

dv r 



This leads to the third condition: 

dC + dC 
dr dv r 



(3.53) 



(3.54) 



The forth condition is not located at the critical point but at the base of the wind, where our solution sho uld fit 
the c o nditio ns of the photosphere. We can specify the optical depth or the wind velocity there. The Eq. 3.39, 



3.52 & 3.54 can be solved analytically for r c , v rc , and v'. CAK solved for M instead of r c . In any case we have 



to assume a value for the forth unknown, find the critical CAK-point, integrate Eq. 3.39, and finally improve 
our assumption about the forth unknown using the boundary condition at the base of the wind. This procedure 
converges sufficiently fast. 

The CAK model gives us a simple description for the momentum transfer from the radiation to the wind. 
But this includes no description for the energy transfer between radiation and wind. Therefore it is not possible 
to solve simultaneously the energy eq uatio n for the local temperature. Instead we have to assume a temperature 
distribution before we can solve Eq. [3.39 for the velocity distribution. A separate model can then be used to 
improve the temperature distribution as CAK did |l6|, This iteration converges fast, because the velocity 
distribution depends only weakly on the temperature distribution. 

When we ask for a wind solution with constant but nonzero temperature, we can apply our discussion about 
the critical CAK-point. So we would expect a unique solution. But if we write down the third critical point 
condition for w 

dA ., dB 



dw 



du 







(3.55) 



dA 
dw 
dB 
du 



2w 2 



> 



= 2vi > 



(3.56) 
(3.57) 



wc sec that w must be negative. This would lead to complex and therefore unphysical solutions for Eq. 3.44 



The basic theory, described here, had been improved in many ways. This includes the improved fit formula 
for the force multiplier by Abbott H , the correction for the finite size of the star Jl9|, |60t , the inclusion of multiple 
photon scattering j38, 23 1, and the detailed line balance models of Kudritzki p2, [59[. But in this thesis we want 
to emphasize the effect of magnetic fields. Therefore we will restrict ourself to this basic formalism. 



Chapter 4 

Magnetic fields in stellar winds 



4.1 Introduction 

Already Parker j5i| argued about the role of the magnetic field in the solar wind. But the solar magnetic field is 
weak. Its contribution to the momentum balance is small compared to the thermal pressure of the hot corona. 
Therefore he neglected the Lorentz force as we did in Chapter ||. But he saw the importance of the magnetic 
field for the angular momentum balance of the sun. 

Weber & Davis (69| were the first who replaced the hydrodynamical treatment of Parker by a magnetohy- 
drodynamic (MHD) treatment to analyze the importance of the magnetic field for the solar wind. The second 
important ingredient for their model was the rotation of the sun, which was neglected by Parker as well. Without 
rotation the material will stream out radially along radial, open field lines. In this case there will be no net force 
and we end up with the physics of Parker's model. 

Biermann & Cassinelli || [| developed an improved version of the Weber & Davis model, which treats the 
magnetic field in the equatorial plane without any further simplifications beyond the basic assumptions we made 
already as well. This allows for an oblique magnetic field at the stellar surface, which occurs in fast magnetic 
rotators. This model is the basis for our analysis of magnetic stellar winds. We present here the derivation 
in a generalized version. Biermann & Cassinelli and Weber & Davis assumed that meridional diameter of the 
fluxsheet is proportional to r. This leads e.g. to radial magnetic field proportional to r~ 2 . We introduce an 
arbitrary function a(r), which allows us to describe the derivation from this r~ 2 law due to the influence of a 
wind compressed disk or a wind, which is diluted in the equatorial plane by a non-equatorial wind collimated 
towards the poles. 



In Sect. 4.2 we derive our wind equations from the fundamental equations of magnetohydrodynamics. In 



Sect. 4.3 we discuss briefly simplified versions of our equations. And finally we discuss in Sect. 4.5 the spin-down 



problem, which gives limitations for the model parameters. 



4.2 The model 

We start our derivation from the magneto-hydrodynamical equations for an ideal, perfectly conducting gas in a 
stationary situation. 

In the case of the Parker wind we had a perfect spherical symmetry, which reduced the problem directly to 
one spacial dimension. The one dimensional solution was valid for all areas in the wind. For a rotating star we 
have only one axis of symmetry. Therefore we can reduce our problem strictly only by one dimension. For a 
stellar wind it is most useful to introduce spherical coordinates parallel to the rotation axis and then to neglect 
the dependence of any physical quantity on the azimuthal angle 4>. For winds from other objects (e.g. disks) 
a cylindrical coordinate system might be more useful. In any case we should solve the MHD equations in the 
two dimensions of the meridional plane. In this plane we can not expect more symmetry because the MHD 
equations mix the radial symmetry of gravitation and radiation with the cylindrical symmetry of the centrifugal 
force. A rigid two dimensional treatment is beyond the scope of this project. But we can reduce the problem 
to one dimension if we assume a priori a certain geometry of the flow in the meridional plane. This concept is 
called the fluxsheet approximation. We have to assume the path of the fluxsheet in the meridional plane and its 
cross section, varying along the fluxsheet. Then we need only to find the various physical quantities along the 
fluxsheet. This is a one dimensional problem. In this chapter we choose the equatorial plane as fluxsheet. From 
observation we know, that our sun has a wind in the equatorial plane. And the principal symmetry between 
northern and southern hemisphere of a star makes this scenario quite probable for other stars as well. The only 
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A) 



B) 





a( r) r 2 








r 









Figure 4.1: An equatorial fluxsheet in the (r, 6) plane. Fig. A) shows the situation in previous models. The 
thickness of the fluxsheet in the meridional plane is proportional to r. Due to the rotational symmetry the same 
is true in the equatorial plane. The cross section of the fluxsheet is therefore proportional to r 2 . Fig B) The 
dimensionless function a(r) describes the deviation of the fluxsheet cross section from r 2 . 



other reasonable configuration in the equatorial plane is an accretion disk. We will discuss non-equatorial winds 
m Chap. 

We start from the equations of magnetohydrodynamics (MHD) for a nonviscous, perfectly conducting, and 
quasi stationary gas: 



V • (pv) 
p{\ ■ V)v 
V x (v x B) 



= 



F ext --Bx(VxB) 
0. 



(4.1) 
(4.2) 
(4.3) 



The external forces 



VP + F r 



ad 



(4.4) 
(4.5) 



include gravity, radiation pressure, and thermal pressure as discussed in the previous chapters. Since we restrict 
our discussion to the equatorial plane {9 — tt/2) we can set ve and Bg to zero. But there derivatives with respect 
to 8 might not be zero. The momentum equation then simplifies to 



dv r 
dr 



Or 

„,2 



h [(V X B) X B] < 



B r drB^ 
Anr dr 



1 



GM 



r 

,2 



r da 
a dr 
1 d 



dvl 
dr 



"5rad H 



r 87rpr 2 dr 

In the fluxsheet concept the conservation laws for magnetic flux and matter are given by 

Fb — B r r 2 a(r) — const 
M = Airpr 2 v r a{r) = const. 



(4.6) 



(4.7) 



(4.8) 
(4.9) 



1G 



CHAPTER 4. MAGNETIC FIELDS IN STELLAR WINDS 



The definition of M can be misleading. It is in analogy to the definition for a spherical wind, where we have 
a(r) = 1. In the case of a cylindrical wind with rotation we must calculate the wind for all latitudes to get an 
accurate value for the total mass loss of the star. But if we assume that the wind is mostly spherical at the base 



of the wind (r = R), and if we normalize a to a(R) = 1, Eq. 4.9 will give a good estimate for the total mass loss 



rate. If we assume that the mass loss is reduced at the poles due to the missing rotation there, Eq. 4.9 will give 



an upper limit for the real mass loss rate. Which case is true, depends on the rotation rate of the star. 



The flux freezing condition Eq. 4.3 for the perfectly conducting wind in the equatorial plane leads to 



d d 
= j- [r (v^B r - VrBj,)} - — (vgB^ - vM . (4.10) 

The symmetry between northern and southern hemisphere requires that all quantities except vg and Bg have a 
vanishing derivative with respect to 9 in the equatorial plane. Now we replace the 0-derivatives of vg and Bg 



using V ■ B = and Eqs. 4.1, 4.8 & 



d o! 

= 7— [r (v^Br - V r Bj,)] (VrBj, - v^B r ) , (4.11) 

or a 

where we used the ' to mark the derivative with respect to r. This equation can be directly integrated and leads 
to 

ar(v r Bj, — Vtf,B r ) — const (4-12) 
= —QciqR 2 B r0 (4.13) 
= -QF B . (4.14) 

Here we introduce the concept of the 'base of the wind'. At the base of the wind we have a negligible radial 
wind velocity. It is the transition region between the static envelope and the dynamically expanding wind of the 
star. We use R for the radius of the base of the wind and denote other physical quantities at the base of the 
wind with the index o- For stars with a weak, optically thin wind, like our sun and most main sequence stars, 
the base of the wind is identical with the photosphere defined by an optical depth of one for electron scattering. 
This is e.g. not true for Wolf-Rayet stars. We denote the rotation rate of the star with fL It is defined by the 
rotation of the quasi stationary magnetic field configuration. At the base of the wind the radial velocity v r o can 
be neglected compared to the co-rotation velocity u ro to = QR- If at the base of the wind the radial component 
B r Q and the azimuthal component B^q of the magnetic field are of the same order, the azimuthal component 
of the wind velocity v^ will be nearly identical to the co-rotation velocity at the base of the wind. In this case 
Q, will be the rotation velocity of the stellar surface as well. We will later find this confirmed in our numerical 
models. It is important to note at this place, that f2 may depend on latitude. This is e.g. the case for our sun. 
From theoretical models we know as well, that the rotation rate can change drastically inside the star. 

Due to the rotational symmetry of our star the angular momentum component parallel to the rotation axis 



is a conserved quantity. Equation 4.6 describes the transport of angular momentum in the wind. Using the 



conservation laws for magnetic flux and matter we can integrate Eq. 4.6 to yield the angular momentum loss of 
the star per unit wind mass 

Br 

L,j = rvA, tBa, — const. (4.15) 

J v 4:TTpV r 

This includes the angular momentum in the wind matter and the angular momentum in the bent magnetic field. 
We parametrize Lj by 

L-j = eilR 2 , (4.16) 

where flR 2 is the angular momentum transport per unit mass by matter at the base of the wind. The angular 
momentum transport of the magnetic field at the base of the wind is included through the dimensionless factor 
e. For all reasonable configurations e is larger than 1. We discuss this later in detail. 
The radial velocity can be expressed relative to the radial Alfvenic velocity by 

V r V r , 

M\ r = — = -B-. (4.17) 



VAr 



This leads to 



Bf a 2 r* _ Fg 1 _ 1 v r0 a R 2 



M\ r Airpv 2 a 2 r 4 M v r ar 2 M\ r0 v r ar 2 
and 

1 r 2 1 F% n 



Ml r R 2 e M L jVr a 

1 v r0 1 a 



M Ar0 



(4.19) 

,> ( 4 - 2 °) 
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With Eq. 4.14,4.15fe4.19 we can now express by 
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where we have introduced the auxiliary quantities 

U e = 1- 

U M = 1- 

= 1 - 



1 



1 ?V0 

Mi r0 « P 



Oo 

a 



Now we have expressed v^ in terms of r and v r . Before we can solve Eq. 4.7 as an ordinary differential 
in v r (r) we have to do the same for 



1 

v r r 



rv<pB r 



VF B 



a 
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v r r 
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nF B 


1 


v r r 
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u t 
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Um_ 
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rB r L, 



riF 
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fir 2 
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Fl Q 1 
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The Lorentz force term in Eq. 4.7 can be split into two parts by 



d 



%ixpr 2 dr 



(rB 4 



av r d 
2M dr ' 

a d fv 2 (rB^) 2 




Each part can now be expressed in terms of r, v r , and v' r . 



nF B u e 



v r M V av r U M , 
Or\ 2 a 1 v r0 R 2 ( U e Y 
v r J a M\ r0 v r r 2 \U M ) 

1 _Oof^R\ 2 fVroY(Ue 



M ArO a \ l 'r0 ) \ V r J \ U M 



(4.21) 

(4.22) 
(4.23) 

(4.24) 
(4.25) 

(4.26) 
equation 

(4.27) 
(4.28) 
(4.29) 
(4.30) 

(4.31) 



(4.32) 



(4.33) 
(4.34) 
(4.35) 



a n 2 F 2 d ( u e 



2Mv r dr \aU M 
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MlrO 



VrO 
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n 2 U 2 
a u M 



(4.36) 
(4.37) 
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(4.38) 



(4.39) 



(4.40) 



Here we used again the ' to mark the derivative with respect to r. It useful to collect terms depending 
because this makes it later easier to find v' as function of r and v r . 



on v' 



I + III = 
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(4.41) 
(4.42) 



As final result for the Lorentz force term we get then 
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This can now be combined with the term for the centrifugal force 
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Now we can write the radial momentum equation (Eq. 4.7) as an ordinary differential equation in v r (r), which 
allows the desired one dimensional solution: 



GM 
1 - 



dv r 
dr 



1 



VrO 



QR 

V r0 



9 

v;r 



GM 
fl 2 R 3 R 1 



GM r U 2 M 



ra 
a 



GM dr 

1 VrQ do 



M 2 Ar0 



T^ffrad 

GM 
'U M 



1 



U^ao 

U M a 



1 VrO dp 



n 2 R 3 r 1 v r0 U 2 a ra' 



' M 



(4.49) 



GM R M 2 rQ v r Ul 

In this chapter we neglect the radiative acceleration although we included the radiation term for reference into 



Eq. 4.49. Here we want to concentrate on the effect of the magnetic field. The CAK model for the line acceleration 



acceleration introduces a radiative term which is nonlinear in v' r . We discuss the consequences of such a term in 
the next chapter. 

For the numerical solution of this equation we introduce dimensionless variables 



V = — 

r* 



GM = 



(4.50) 

(4.51) 
(4.52) 



Using u ~ r^ 1 instead of u ~ r has the advantage that the region close to the stellar surface, where the wind 
is strongly accelerated on a short length scale is enlarged. This region is the interesting region, where most 
things happen, the enlargement in the u-scale will ensure proper solutions for the wind equation in this region. 
On the other side we have the region far away from the star, where the wind is streaming outward with a 
nearly constant velocity. This region is physically and numerically uninteresting. In the u-scale this region is 
shrunk. And therefore not much numerical effort will be spent to solve the wind equation there. With a proper 
temperature stratification this region can even extend to r = oo (u = 0). Chap. ^ shows several plots with wind 
solutions as functions of u, which show nicely that u is the appropriate coordinate to compute and to plot wind 
solutions up to large radii. Our choice for the length and velocity scale has still one degree of freedom. We can 
still choose an arbitrary value for r* or v*. This of course has no influence on the physics. But an improper 
choice can trouble the numerics. One good choice for fast magnetic rotator models (strong magnetic field, fast 
rotation) is the Michel velocity 
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We choose for simplicity just the stellar radius as length scale 



This length scale is reasonable for all wind models 
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= R. (4.54) 
even for nonrotating or nonmagnetic wind models. The 



dimensionless version of Eq. 4.49 
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where we have used the dimensionless quantities: 
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In analogy to Eq. 3.39 we can write Eq. 1.55] as 

A(u,y)y + B(u,y) = 0. 



(4.57) 
(4.58) 
(4.59) 
(4.60) 
(4.61) 
(4.62) 

(4.63) 



Equation 4.55 is the fluxsheet generalization of the results of Bicrmann & Cassinelli M and Bicrmann 



4.3 Simplified versions 

For this model it is important that both, magnetic field and rotation, are present. If we would turn the rotation 
off (A = ?/ r oto = 0) Eq. 4.49| directly simplifies to the equation of the Parker wind. The open magnetic field lines 
would be stretched out radially by the wind. 

If we take rotation into account but assume that the star has no relevant magnetic field, the conservation of 



angular momentum (Eq. 4.15) would simplify to 

r"U0 = const 



For Eq. 4.55 we find then 
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(4.65) 



This scenario in combination with a CAK radiation force has been analyzed by Pauldrach et al. 



4.4 The critical points 



As in the case of the Parker or CAK wind it is important to understand the critical points of Eq. 4.55 to find 



proper physical solutions. In this model we have not one but three critical points. Under certain circumstances 
the solution need not pass through all these points. But nevertheless it is important to understand them all. 

The most important critical point is the Alfvenic point. At this point with r = ta c the radial wind velocity 
equals the radial Alfvenic velocity (M^ r — 1). And therefore Um will be zero at this point and change its sign. 
If we require that B<f, is finite and continuous, U e must be zero and change its sign there as well (Eq. |4.3C| ). For 
the radius of the Alfvenic point and the radial velocity there we find 
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As in Chap. @ we can try to use l'Hospital's rule again to find y at the critical point: 
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(4.66) 
(4.67) 

(4.68) 



In the case of the Weber & Davis wind we do not find a finite number of possible values for j/ac- Therefore we 
do not have a simple X-type critical point here like in the Parker wind. This local analysis allows a range of 
values for ijAc- We find in general 

Ue _ Bff, V r 

Um B r fir 



(4.69) 
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Depending of the values of v r and B, U £ /Um can have different values. This leads to different values for y Ac 



according to Eq. 4.68 



Since U c /Um is finite, and Um equals 0, at r Ac the U 2 /U^ term in A(u, y) will dominate for r ~ r Ac - For 
r ^ r Ac we have A > 0. But if the velocity at the base of the wind is small enough, A(uo,yo) will be less 
than zero due to the then dominating thermal pressure term. For r > r Ac we have A < 0. And for sufficient 
high velocities at r = oo A will be positive due to the inertia term. So we have two radii R < r s < r Ac and 
f"Ac < < oo where A will be zero. At these radii we have a X-type singularity as in the Parker wind. We 
denote these critical points and their analogies in the case of an magnetic CAK wind as the inner and the outer 
critical points. The condition A — leads to 
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where we have used the total Alfvenic velocity 



Equation 4.72 is the well known dispersion relation for slow and fast magnetosonic waves |27j. Therefore at the 
slow and fast critical point the radial wind velocity v r equals the local phase velocity for slow respectively fast 
magnetosonic waves. There is a deep physical connection to the problem of causality in a stationary physical 



model. Weber & Davis |69|, Figs. 1&2] described the topology for the solution Eq. 4.63 and the role of the critical 
points. In the next chapter we combine the Weber & Davis model in the form presented here with the CAK 
model. There we will discuss the critical points and the solution topology in detail. 



4.5 The spin-down problem 

One important aspect in the theory of magnetic winds is the angular momentum balance. Every star which 
rotates and loses mass will lose angular momentum as well, because the rotating surface layer of the star carries 
a certain fraction of the stellar angular momentum. When this layer evaporates into the wind it will take its 
angular momentum along. Due to a magnetic field or viscosity the angular momentum loss may be enhanced. 
We expressed this by the factor e in our formula for the angular momentum loss per unit mass (Eq. 4.16| ). A 
normal star without an enhanced angular momentum loss has e = 1 . The angular momentum loss rate^ of a star 
is 

dM 



L = I Li-fa**' ( 475 ) 



where Lj is the angular momentum loss per unit mass (Eq. 4.15) and u> is the solid angle. In order to evaluate 
Eq. 4.75| properly we would have to calculate a whole 3-dimcnsional wind model. This is is beyond this thesis. 



So we have to assume the star is sufficiently spherically symmetric at its surface, so that we can use the results 
from our equatorial wind model. This assumption is only good if the rotation rate is sufficiently low. Since 
this is not the case for most of our models, the results of this section should only be used carefully as a rough 



estimate. We find under these assumptions from Eq. 4.75 



L = -eMR 2 n = -Mr 2 Ac n. (4.76) 
3 3 

The interesting point now is to compare the timescales for mass loss 

r M = | (4.77) 

and angular momentum loss 

tl = |. (4.78) 

1 L = dL/dt is beside M the second exception from the rule that the dot ' marks a derivative with respect to u. 
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B-field 



Lorentz force 




B-field 




B-field inside the star 



A) 8 > 1 



B) 8 < 1 



Figure 4.2: 
e < 1 (B). 



shows the magnetic field configuration in the equatorial plane of wind models with e > 1 (A) and 



We know from observation and stellar evolution models that massive stars lose a large fraction of their mass 
before they explode as a supernova. I.e. their mass loss timescale is comparable to their lifetime. If rotation 
plays an important role through the whole life of a massive star, its angular momentum loss timescale should 
not be much smaller. Otherwise the star will spin-down already in an early stage of its evolution. The problem 
now is to estimate L, the total angular momentum of our star. Today no models for the inner structure of 
rapidly rotating stars with significant magnetic fields are available. The most simple by extremely unphysical 
assumption is a rigidly (!) rotating star of constant (!) density. In this case we have L = {5/2)MR 2 n and 
therefore 

5 e 

This result suggests that every star has an angular momentum loss timescale smaller than the mass loss timescale. 
The naive interpretation would be that even a star without magnetic field (e = 1) would spin down rapidly. But 
of course, a star without angular momentum transport between different layers of material will not spin down 
at all, if we assume that it does not expand or shrink. If the angular momentum is conserved in every layer 
then every layer will keep to overall rotation rate. Due to all the simplifications Eq. 4.79 allows to criticize wind 
models only if tm and tl differ by more than one or two orders of magnitude. A final answer for the spin-down 
problem requires a unified model which describes the evolution of the star and its wind for all latitudes in detail 
over the whole lifetime. This is beyond todays possibilities. 

Therefore our wind models are reasonable as long as they have a small value of e. But what is the physical 
meaning of e? It is the factor by which the angular momentum loss is enhanced compared with the case of no 
magnetic field. From our discussion in the previous section we see that the radius of the Alfvenic point is given 
by r^ c = eR 2 . If we make the quite physical assumption that the magnetic field transports angular momentum 
from the star away, we will always have e > 1 and therefore r^ c > R. A star with e less than 1 would lose less 
angular momentum than a star without magnetic field. If we assume that at the stellar surface (or somewhere 
just below) we have a layer of stellar material without significant outward motion (v r = 0), we will have strict 
co-rotation there due to the flux freezing. When this layer evaporates into the wind, the magnetic field would 
have to transport angular momentum back into the star to keep the angular momentum loss less than in the case 
of a nonmagnetic star. Since the nonmagnetic star rotates with a constant rotation rate, a magnetic star with e 
less than 1 would spin up like a rotating firework. Figure [4.2| shows shows the magnetic field configuration in the 
equatorial plane of wind models with e > 1 (A) and e < 1 (B). In the first case (A) the magnetic field line pushes 
the wind material outside and in the direction of the rotation. The first effect enhances the wind outflow. The 
second effect transfers extra angular momentum from the magnetic field to the wind matter. This extra angular 
momentum comes via the curvature of the field line from the star. In the second case the magnetic field line is 
mostly tangential at the stellar surface. The wind accelerated by radiation pressure pushes against the field line. 
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This bends the field line clockwise. Through this curvature angular momentum is transferred from the wind 
material via the field line back to the star. Therefore we have e less than 1. The consequence of the curved field 
line is a Lorentz force which decelerates the wind in radial and azimuthal direction. If this effect is very strong 
the wind can even start to rotate counterclockwise. In any case scenarios of type B have a decelerating Lorentz 
force close to the stellar surface. The force balance close to the surface fixes the mass loss rate. Therefore wind 
models of type B will have a lower mass loss rate than the corresponding model without magnetic field. It is 
not even clear whether a configuration of type B is stable, because the field line has to be bent clockwise inside 
the star (dotted line in Fig. 4.2.B). Otherwise it would hit the clockwise bent stellar surface again. This would 
violate our assumption of rotational symmetry. We therefore conclude that e must be larger than 1 in physical 
wind models with magnetic fields. And we will restrict our further analysis to that case. 



Chapter 5 

Magnetic CAK-winds 



5.1 Introduction 



In this chapter we discuss various numerical solutions for the equatorial wind model developed in the previous 
chapters. Where possible we compare our results with the literature. The combined wind equation for line driven 
winds with magnetic fields is 

A(u, y)y + B(u, y) - C(«, y)y a ™\ (5.1) 

We know this form of the wind equation already from our discussion of nonmagnetic, line driven winds in Chap. [3] 
(Eq. p.3S| ). But now we use the dimensionless quantities (Eqs. 4.5C-4.52) and extend the coefficients A and B 
by the magnetic terms derived in Chap. || (Eq. [4.55 ): 
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(5.3) 
(5.4) 



Eq. 5.4 distinguishes from Eq. 3.43 beside the overall factor r 2 / G M by t he geom etry factor a(u) acak . This factor 
comes from the equation of continuity, which we used between Eq. |3.42| and Eq. 3.43. The geometry factor a(u) 
was not introduced before chapter 0. 



A solution of this wind equation can pass through up to three critical points. In Sec. 4.4 we discussed the 
three critical points of a magnetic wind without line acceleration (i.e. C = 0). The Alfvenic critical point, we 
discussed there, exists in a line driven, magnetic wind as well. And our arguments, that A will pass through zero 
somewhere inside and outside the Alfvenic point, are still valid. But since the we have now the line driving term 
Cy a ° ak in our Eq. 5.1, we can not use the old Parker argument that A(u,y) = implies B(u,y) = 0. Therefore 
we do not have the classical X-type critical point at the radii where A(u, y) is zero. Rather A and B change 
their signs in the subtle way we described for the critical point of nonmagnetic CAK-wind in Sec. 3.4. So we 
have two critical points of the CAK-type. 

In order to compute a numerical model we have to fix several model parameters. The star is described by its 
mass M, its luminosity L, its radius R, its rotation rate Q, and its surface radial magnetic field strength B r $. 
For the wind we have to specify the CAK-parameters a ca k and fc C ak, the temperature profile T(r), the electron 
scattering opacity KTh, and the meridional shape of the wind described here by the area function a(r). When 
we put all these numbers into Eqs. 5.1-5.4 we find that two further quantities remain to be fixed. E.g. we can 
specify additionally the mass loss rate M and the radius of the Alfvenic point tac- These two parameters are 
very crucial because t hey a llow to relate our numerical models to observation and to control the stellar angular 
momentum loss (Eq. 4.76 ), which is important for the stellar evolution. Therefore it would be nice, if these 
parameters were rather results than input parameters of our models. This can be realized, if we require that 
our wind solution has to pass through all three critical points. Through this additional assumption M and r& c 
become eigenvalues of our wind equation. We find these eigenvalues using the shooting method. For given values 
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of M and ta c we can find the two CAK- typ e points by solving twice the nonl i near set of equations . 3.52| fc 3.54 



From these poi nts w e can integrate Eq. 5.1 to the Alfvenic point (Eqs. 4.66 ^: 4.67 ), to the stellar surface, and to 
infinity. In Sec. 4.4 we argued that the acceleration y is not fixed a priori at the Alfvenic critical point. But, of 
course, y should be continuous at the Alfvenic point. Therefore we should find the same value of j/(uac) when we 
integrate Eq. 5.1 from the inner and from the outer CAK-point to the Alfvenic point. This is the first condition 
for the eigenvalues M and ta c - The second condition is related to the stellar surface. We defined a stellar radius 
R as the radius at the base of the wind. In our models we require a certain, small, initial wind velocity yo at this 
radius. This condition is very reasonable from the theoretical point of view, because a theoretical wind model 
should describe, how the wind material is accelerated from the low, subsonic velocity of the quasi static, stellar 
atmosphere to the high observed terminal velocities. Alternatively we could require a certain optical depth (e.g. 
To = 1) at the stellar radius R. For stars with an optical thin wind like O and B stars these two alternative 
definitions of R are practical identical. For stars with an optical thick wind like Wolf-Rayet stars the optical 
depth at the subsonic base of the wind is not known. If we would start our wind model for a Wolf-Rayet star at 
r = 1, we would describe only the outer, high velocity part of the wind. Such a model could not be considered as 
a complete wind model. Additionally specifying r requires more numerical effort than specifying y a . Therefore 
we prefer t o sp ecify yo. 

In Sec. |5.2| we compare our res ults with the model of Friend & MacGregor (20) , who published a model with 
a different wind equation. In Sec. 5J3 we our model to the case of only two critical points. This is important for 
our study of winds outside the equatorial plane in Chap. [?]. 



5.2 The wind models of Friend &; MacGregor revisited 

We start our discussion with wind models which have no compression or expansion of the wind in meridional 



direction, i.e. the cross section function a(r) as defined in Eqs. 4.8&T.E is constant. Additionally we assume 
that the wind passes through all three critical points. Such a model can be compared to the model published 
by Friend & MacGregor p0[ . Their model is based on the CAK model for the radiative acceleration and the 
Weber Davis model for the magnetic acceleration (their Eqs. 1-21). Up to this point their model is identical 
with our model. But from these equations they derive a final differential equation for v r (their Eqs. 22-26) which 



differs in the magnetic terms from our final equation (Eqs. 5.1-5.4). And this is not only a question of different 
notation. Our final equation is derived strictly from the common basic assumptions (their Eqs. 1-21) without 
any further assumptions or approximations. Due to the necessary brevity in an article Friend & MacGregor 
give no further information how they obtain their final equation from the basic equations. Therefore we can 
not analyze the differences between the final equation of Friend & MacGregor and our final equation in further 
details. To study the differences of our models in the numerical level we recalculated their numerical models 
with our software using our final wind equation. They calculated 16 models with varying magnetic field and 
rotation rate for the 06ef star A Cephei. For this star they used the parameters: M = 50M Q , L = 6.76 x 10 5 Lq, 
R = 19.7i?0. The sound speed v s = 30km/s is constant in the wind. This leads to a wind temperature of 
about 65000 K. For the radiative acceleration they claim to use the original CAK Jl6l values a ca k = 0.7 and 



fceak — 1/30. From the statement that the CAK terminal velocity (their Eq. 30, our Eq. |3~5ll) is 1240 km/s we 



can derive a value of T = 0.319 for the Eddington factor. Therefore we find for the electron scattering opacity 
per unit mass Kih = 0.309 cm 2 /g (called a e by Friend & MacGregor). For the CAK mass loss rate (their Eq. 31, 



our Eq. 3.50 ) they quote a value of 5.2 x 10~ 6 M Q /yr. For their parameters we find a value of 9.74 x lO _7 M0/yr. 
We find the same discrepancy in the mass loss rates of the numerical models. We can fix this by assuming a 
modified CAK parameter fc ca k = 0.107 instead of fc ca k = 1/30. Finally we can deduce from their Tab. 2.B that 
their reference radius for the terminal velocity of their numerical models is 200 times the Alfvenic radius. For 
all numerical models we state the rotation rate in units of the critical rotation rate 

u rot Q QR c s 

«rot = = = , (5-5) 

«rot,crit "crit ./<2M(X — T) 

ctrot gives a good intuitive feeling, whether the star rotates fast or slow. In the case a ro t = 1 the centrifugal and 
radiative (Thompson scattering) forces will compensate gravity at the equator. An eruptive mass loss would be 
the consequence. This case is called the J7-limit by Langer & Heger |34| ]. They claim that stars can really reach 
this limit. 



In Table |5.l| we display analogously to Friend & MacGregor's Tab. 1 our results for the 16 wind models 
extended by our values for the wind efficiency and the ratio between terminal velocity and the total Alfvenic 
velocity. We find larger differences only for the radius of the slow (« 7%) and fast (« 20%) critical points. In 
the first case this is due to our slightly different inner boundary condition. We require a fixed radial velocity 
(iVo = v s /20) at the photosphere. Friend & MacGregor required an optical depth for electron scattering of 
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Table 5.1: Wind models for A Cephei analogous to the models of Friend & MacGregor pOj . 
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t cs = 1 at the base of the wind. In the second case th e di fferences are due to the bad numerical conditioning for 
finding the far critical point. All other values in Tab. 5T differ less than 3% between Friend & MacGregor's and 
our calculations. This difference is too small to distinguish any physical from numerical effects. Therefore we 
have to conclude that we found the same wind solutions as Friend & MacGregor although we used a different 
wind equation. Their additional assumptions or simplifications, which lead to their Eqs. 22-26, appear to be 
well chosen. For a deeper physical interpretation of these results we refer the reader to the paper of Friend 
& MacGregor, since we do not want to reprint all their arguments here in detail. In order to check our inner 
boundary condition we calculated models with different values of v r o up to 3v s . We found that this has only a 
very small influence on the wind models. The model parameter like M and changed not more than 10%. This 
shows that our boundary condition at the base of the wind is not critical for our results. For most models they 
changed only a few percent. This weak dependence on the inner boundary condition explains also why we got 
the same results as Friend and MacGregor although they use a completely different inner boundary condition. 
Figs. |5 . l| — [5 . 4j show the velocity profiles for four of our models plotted as function of R/r (the white curve). These 



figures show the different regions in the (u = R/r,v r ) plane, which have zero, one, or two solutions for Eq. 5.1. 
The position of the three critical points are marked by white points. Tracing the solution from the stellar surface 
(R/r = 1) to infinity (R/r = 0) we can see how A(u,y) and B(u,y) evolve through the wind and why the 
solution has to pass through the two CAK-type critical points, if we require, that the wind starts at a subsonic 
velocity and extends to arbitrary large radii. 

The interesting point in this section is that the approximations made by Friend & MacGregor leading to 
their version of the wind equation seems to be quite appropriate. It is not clear whether this is still true in the 
case of higher rotation rates. We will use higher rotation rates in Chap. for our non-equatorial models. 



5.3 Models with only two critical points 



5.1 



1.4 



has two solutions extends down to the stellar 
we argued that the wind has to pass through 



When we look at Fig. 5.5 we see that the region where Eq 
surface - but only with a supersonic wind velocity. In Sec. 
the inner critical point, because only if we switch the solution branch at the inner critical point the solution can 
extend into the subsonic part. But it is not clear whether our wind equation is still valid in the subsonic regime. 
Effects like turbulence, optical thickness, pulsation, meridional motions or other might dominate the wind there. 
In this case we might get wrong results if we try to describe the wind down to the subsonic regime with our 
model. Nevertheless if we restrict our calculations to the supersonic regime we still have to start close to the 
star with a wind velocity small compared to the terminal velocity. Otherwise our model is irrelevant because it 
does not answer the fundamental questions of wind physics: (1) How does the wind material escape from the 
gravity of the star? And (2) how is it accelerated to the high observed terminal velocities? 

To analyze this question we modified our program, so that the wind solutions pass only through the outer 
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0.4 0.6 

u — R/r 

Figure 5.1: shows the radial wind velocity v r (white curve) as function of u = R/r for our version of Friend & 
MacGregor's model 1 (weak B field & slow rotation). The three critical points (from left to right: the outer, 
the Alfvenic, and the inner critical point) are marked by white dots. Additionally we plot, using different gray 
scales, for every point in the (u,v r ) plane the local solution logic for Eq. 5.1, which is the magnetic version of 
Eq. 3.39. The "local solution logic" is explained for the nonmagnetic case in Sect. EO| and Fig. 3.1. The different 



cases in Fig. 3.1 are connected to the different grayscales in this plot by the legend on the right hand side. In the 



magnetic case we have two crit ical points of the (nonmagnetic) CAK-type: The inner and outer critical points. 
At the Alfvenic point Fig. 3T.F applies. The wind solution shown in thi s plot passes from the stellar surface 
(u = 1) to infinity (u = 0) through the following sequence (c.f. Fig 3.1: C-E-A-B(here is the inner critical 
point)-A-D-F(here is the Alfvenic critical point)-C-A-B(here is the outer critical point)-A-D. The situation 
around the Alfvenic point can be seen better in Fig. ETa. 



and the Alfvenic critical points. Through this step we lose one of our conditions for the eigenvalues of Eq. 5.1. 
Now we have to specif y ei ther the acceleration at the Alfvenic point or the wind velocity at the stellar surface in 
order to integrate Eq. 5.1 between the stellar surface and the Alfvenic point. Therefore we do not find a unique 
solution any more. Rather we find a set of solutions with M as a function of ta c or vice versa. 

We calculated solutions for the four corner models of Friend & MacGregor (models 1,4,13,16). Hereby we 
used an initial wind velocity at the stellar surface of v r o = 1.3 v s — 39 km/s. In Tab. |5.2| we list for the four 
models the solutions with the smallest and the largest radius of the Alfvenic point (rAc)- As reference level for 
the terminal velocity we chose 100i?. This allows us to compare these results with the models in the next 
paragraph. For the models 1,4&13 we find that the Alfvenic radius ta c is reduced compared to the models 
with 3 critical points. But for model 16 we find only solutions with an enlarged Alfvenic radius. All other 
quantities vary up to an factor of two. This shows that even for fixed magnetic field strength and rotation rate 
magnetic wind models have a strong capability to fit observational constraints. Additionally we see that a proper 
treatment of the wind physics at the stellar surface is crucial for quantitative wind models. We did not expect 
this in the last section where we found that models with 3 critical points are insensitive for the in itial wind 
velocity u r o- Therefore it is important to improve the wind physics near to the stellar surface. Figure 5.7 shows 
our wind solution 13a. The interesting aspect about solution 13a is that it has a significantly reduced angular 
momentum loss (e = 2) compared to the model with 3 critical points no. 13 (e = 3). Additionally model 13a 
has a much steeper velocity law and a much higher terminal velocity than model 13. This can also help to fit 
observations. In Fig. 5.7 we show the velocity profile for our wind model 13a. 

It is also possible to question whether the outer critical point exists. It will certainly exists when we have 
a finite temperature at an sufficiently large radius. But is is not necessary that our wind solution extends to 
infin ity in the mathematical sense. This would not be physical as we discussed in Chap. |^. We can see in 
Fig. 5.6 that the wind solution has to pass through the outer critical point before it reaches at r as 200i? a 



region, where again only local solution of Eq. 5T exists. This was the argument which required the outer critical 
point. But we will discuss in the next chapter a model for linear waves in magnetic stellar winds. We will find 
there, that these waves are highly instable and will grow to linear shocks rapidly. Although the latter effect 
demands for a numerical model for nonlinear waves, we can seriously expect that the wind at large radii will 
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Figure 5.3: shows the same as Fig. 5.1 but for our version of Friend & MacGregor's model 13 (weak B field & 
fast rotation). For the technical explanation of this plot check Fig. 5.1. 
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Figure 5.5: shows the inner part of our wind model 16 (Fig. 5.4). The region with two solutions extends down 
to the stellar surface at u = 1. The wind solution can only extend into the subsoni c pa rt, after it has passed 
through the inner critical point. For the technical explanation of this plot check Fig. 5.1 . 
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Figure 5.6: shows the outer part of our wind model 16 (Fig. 5.4). The wind solution can only extend beyond 
u w 0.05 (r ~ 200i?) after it has passed though the outer critical point. For the technical explanation of this 



plot check Fig. 5.1 



Table 5.2: Wind models for A Cephei which do not pass through the inner critical point 



no. 


(G) 


Qirot 


r C 2 
R 


R 


V Ac 

(km/s) 


(km/s) 


10 6 M 
(M /yr) 


M«oo 

L/c 


Hop 
VAoa 


la 


200 


0.218 


1.60 


1.54 
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1751 


4.92 


0.624 


11.89 


lb 


200 


0.218 


2.26 


1.67 


820 


1379 


5.24 


0.520 


8.58 


4a 


1600 


0.218 


9.09 


8.16 


4080 


4384 


2.81 


0.894 


4.45 


4b 


1600 


0.218 


24.26 


9.43 


1684 


1977 


5.09 


0.720 


1.81 


13a 


200 


0.697 


1.50 


1.41 


1611 


2985 


3.73 


0.814 


7.28 


13b 


200 


0.697 


5.04 


1.64 


753 


1380 


5.90 


0.594 


2.91 


16a 


1600 


0.697 


49.05 


7.09 


2220 


3030 


6.84 


1.513 


1.30 


16b 


1600 


0.697 


25.66 


7.84 


3685 


4367 


3.37 


1.075 


1.56 



not be a continuous, smooth, stationary flow, as assumed in the wind model used in this chapter, but rather a 
chaotic system of shocks and nonlinear waves. Therefore it might be that the smooth wind model used in this 
chapter is not a proper description for the stellar wind at large radii - even for the time averaged properties of 
the wind. If we restrict our wind model to smaller radii with r < 200i?, the need for the outer critical point is 
relieved. Therefore we calculate now a set of wind models which have an inner but no outer critical point. As in 
the previous case we have therefore one free parameter in the model. Analogous to Tab. |5.2| Tab. shows our 
results for wind models without an outer critical point. Hereby we restricted ourself to a maximum of 15i? for 
the radius of the Alfvenic critical point r^ c - Especially for wind model 13 such a high value for rj± c is already 



unphysical. Figure 5.8 shows the extreme case of solution 13c. The wind solution starts at the stellar surface 
well above the region of no wind solution. Since our solutions pass through the inner critical point, we could 
have chosen a subsonic wind velocity a the stellar surface. But in order to compare these solutions with our 
solutions from the last paragraph we chose the same initial wind velocity of v r o = 39 km/s. Our solutions pass 
through the inner and the Alfvenic critical point, and stops at 100 stellar radii. This is well below the region 
where the passage through the outer critical point becomes important (r » 300i?). 

Since a picture can tell m ore t han 1000 words we show finally for our wind model 13 the mass loss rate M vs. 
the terminal v eloci ty Woo (Fig. |5.9| ) and the wind efficiency (Mvoo)/ (L/c) vs. the angular momentum enhancement 
factor e (Fig. 5.10| ). The wind model with three critical points seems to form a situation between the two cases 
with only two critical points. From the view point of avoiding the spin-down problem in a strong wind models 
without an inner critical point are very interesting. They allow an increased wind efficiency together with a 
reduced angular momentum loss. This was made possible by the fact that we modified our model close to the 
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Figure 5.7: shows our wind model 13b. This solution has no inner critical point. Therefore the sequence of the 



solution logic from the stellar surface (it = 1) to infinity (u — 0) is (c.f. Fig. 3.1): A-D-F(here is the Alfvenic 
critical point)-C-A-B(here is the outer critical point)-A-D. For the technical explanation of this plot check 
Fig. [O. 





Figure 5.8: shows our wind model 13c. This solution has no outer critical p oint . Therefore the sequence of the 
solution logic from the stellar surface (u = 1) to infinity (u = 0) is (c.f. Fig. 3.1): A-B(hcrc is the inner cri tica l 
point)-A-D-F(here is the Alfvenic critical point)-C-A. For the technical explanation of this plot check Fig. 5.1 



Table 5.3: Wind models for A Cephei which do not pass through the outer critical point 
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1473 
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4d 
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0.250 


0.39 
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1.05 
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0.385 
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0.004 


0.02 
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1.07 


15.00 
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648 


5.27 


0.249 


0.11 
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star. Our modification to assume that the wind is already slightly supersonic at the base of the wind is arbitrary. 
Somehow the wind has to become supersonic. This problem could be solved by using a different theory for the 
wind close to the stellar surface. Such a theory could e.g. incorporate the consequences of an optically thick 
wind. Further out we should smoothly switch to our wind model taking advantage of the high wind efficiency 
and the low angular momentum loss of our wind solutions without an inner critical point. Therefore such an 
alternative model must include a consistent transition from the inner to the outer regions of the wind. 

We did not make yet use of the capability of our theory to describe winds which are compressed or diluted 
in the equatorial plane. We will discuss this aspect in a more generalized framework in Chap. pL 
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Figure 5.10: shows the wind efficiency (Mw (X) )/(L/c) vs. the angular momentum enhancement lactor 



Chapter 6 

Waves and shocks in stellar winds 



This chapter is based on an article 67 by the author and his academic supervisor Prof. P. L. Biermann. 



6.1 Introduction 

In the previous chapters we asked for a stationary solution of the wind problem. This assumption leads to a 
strong simplification of the model. Such a simplified model is a good starting point for understanding a physical 
system by describing its time averaged features. Therefore the first models which were applied to the solar wind 
were stationary (c.f. Chaps. as well. But from observations we know today that the solar wind is far from 

being stationary. The strong fluctuations in the solar corona heat the corona to the high temperatures, which 
are necessary to drive the solar wind by thermal pressure. Therefore even the time averaged features of the solar 
wind can not be understood without taking fluctuations into account. Therefore it is quite possible that the 
wind of hot stars can not be described properly without taking fluctuations into account. 

The analysis of wind fluctuations breaks into three parts: (1) How are fluctuations generated? And how do 
they develop in the stellar wind? (2) Which consequences do they have for the overall model of the wind? And 
(3) which influence docs fluctuations have on observations? 

Wind fluctuations can have two consequences for observation: (1) Is it possible to observe fluctuations directly 
by analyzing the time dependence of observational parameters like line profiles? (2) How do fluctuations influence 
the derivation of physical parameters like mass loss rate and terminal velocity from observational parameters 
like spectra and fluxes? Since this thesis is a theoretical orientated work, we do not want to go to deep into 
observational details here. Observed fluctuations can have two reasons: (1) Periodic fluctuations are probably 
connected to the stellar rotation. Some localized feature on the rotating stellar surface perturbes the wind in the 
line of sight periodically. This scenario is similar to a pulsar. Such fluctuations have e.g. be observed by Gagne 
et al. |21jj . These observations allow to measure the stellar rotation rate, which is crucial for our models. Gagne 
et al. call 9 1 Orionis an 'oblique magnetic rotator' with a rotation period of 15.422 days. This star therefore 
rotates with 3% of its critical rotation rate. (2) The second type of fluctuations have a stochastic characteristic. 
This type of fluctuations are presumably created in the wind by some kind of plasma instability. We will discuss 
such instabilities in this chapter. Further evidence for fluctuations are the observed X-ray and nonthermal radio 
emission. Lucy |}6| proposed that the observed X-rays are emitted by very hot plasma, which is heated up by 
shocks. Shocks are required as well to explain the nonthermal radio emission observed by Bieging et al. || from 
OB stars and by Abbott et al. Q from Wolf-Rayet stars. This nonthermal radio emission is in fact synchrotron 
emission by relativistic electrons. Therefore this gives clear evidence for non-negligible magnetic fields on these 
stars. The relativistic electrons can be produced by Fermi acceleration in shocks. In binary systems these shocks 
could be created by the colliding winds of the two components. For single stars intrinsic wind instabilities are 
the only possible source for the required shocks. 

Fluctuations in the wind influence the observation of basic wind parameters like the terminal velocity and 
the mass loss rate. The observed terminal velocity is inferred from the blue edge of P-Cygni lines. This is in 
fact the maximal velocity reached by a significant fraction of the wind material. In the case of an permanently 
accelerating stationary wind this is approximately the theoretical terminal velocity Woo = v(r — oo). But if 
we have fluctuations in the wind where regions of temporarily enhanced wind velocity have simultaneously an 
enhanced density, we will see the enhanced velocity as the blue edge of the P-Cygni lines. In the linear theory 
such simultaneously enhancements are described by outward running waves. This enhanced velocity is in the 
case of linear waves the unperturbed plus the velocity amplitude of the waves or, if outward running waves 
have already steepened into outward running shocks, the unperturbed Voo roughly plus the speed of the shock 
front. The latter scenario was used by Lucy to explain the X-ray emission from these stars. 
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The observed values for the mass loss rate M are influenced by waves in two ways: (i) M is inferred from 
observed values of the density p by M ~ pr 2 Voa. An overestimated Voq therefore leads to an overestimated M . 
(ii) The values for p are inferred from observed radio or UV fluxes. A matter distribution perturbed by waves or 
shocks generates higher fluxes, which lead to higher values for p if the perturbations are not taken into account 
properly (Abbott et al. || and Hillier pi|). The nonthermal contribution to the radio emission additionally 
increases the estimates for p, if it is not properly separated from the thermal emission. Fluctuations can also 
be an additional source and transport mechanism for momentum. When waves are amplified by radiation they 
extract extra momentum from the radiation and therefore increase the overall efficiency of radiation pressure. 
When these linear waves later dissipate or steepen into shocks, they transfer their momentum to the wind. This 
effect was analyzed by Koninx pfj| . Koninx model and our discussion about the derivation of observed values 
for M and Woo are valid only if the wind is dominated by outward running waves. 

The fluctuations in the solar wind are generated in the convective envelope of the sun. From there they 
propagate into the wind. At the base of the wind they heat up the plasma to the very high temperatures 
(~ f0 6 K) of the corona. Hot stars do not have a convective envelope. Their envelope is mostly radiative. Such 
a radiative envelope will not produce strong fluctuations and therefore presumably no corona. Therefore we 
need another mechanism to produce fluctuations. A good candidate is the radiative line driving force. This 
force is the dominating driving mechanism for the unperturbed wind. Any instability in this force could create 
significant fluctuations. The possibility for such an instability was already mentioned by Castor, Abbott and 



Klein in their famous paper about line driven winds |16 



This idea was later elaborated by several authors. MacGregor et al. |4fJ developed a model for short, optical 
thin waves. They found that these waves are unstable and grow rapidly. This result was opposed by Abbott (!]], 
who found stable waves assuming long wavelengths. Owocki and Rybicki j5^| unified these contradicting models 
by a more elaborate, but still analytical linear model for all wavelengths. This model showed that waves with 
short wavelength (A < Ls) will grow rapidly. They refined their model later in a series of papers ]53|, |54], [S3, 55 1. 
The basic results is that in the thin part of the wind fluctuations are strongly amplified. The fluctuations are 
dominated by inward running waves steepening into shells of increased velocity but decreased density. These 
fluctuations are advected outward by the average wind. These papers all do a linear analysis of the fluctuations. 
Krolik & Raymond |3l| and MacFarlane & Cassinelli |4(| analyzed the structure of a single outward running 
shock shell in a nonmagnetic wind and found that the shock shell is driven by the radiation field to a velocity 
much higher than the velocity of sound and therefore much higher than the phase velocity for outward running 
linear waves found by Owocki & Rybicki |fi^| . But later Owocki et al. [f56| found in a nonlinear time dependent 
calculation, that waves in the stellar wind are dominantly running inward and steepen into inward running, 
reverse shocks. This confirmed their previous result from linear analysis. 

The aim of this chapter is to show in the limit of a linear analysis that the objections of Owocki et al. J5(| 
do not apply, if a strong magnetic field and rotation are present. Therefore we connect the ideas of a rotating 
magnetic field, developed in the previous chapters, and of unstable waves (Owocki & Rybicki |52| in the wind of 
a massive star and do a linear stability analysis. We show that the magnetic field increases the number of wave 
modes and changes their properties significantly. High phase velocities can be achieved due to the high Alfven 
velocity. Therefore inward running waves will not be advected outward and will not steepen into reverse shocks 
anymore, which dominate at large radii. Furthermore outward and inward running waves have the same growth 
timescale in the short wavelength regime where both modes are unstable. So we can expect outward running 
waves far from the star and forward shocks in the nonlinear regime. 

In Sect. 6/2 we outline the model of Owocki & Rybicki [[^J for the radiative instability. In Sect. 6.3 we derive 
the dispersion relation for radiatively amplified waves in the presence of a magnetic field. In Sect. 3.4 we discuss 



our unperturbed wind models. In Sect. p.5| we discuss the waves we found for the wind models of Sect. 6.2. In 



Sect. 6.6 we discuss the observational consequences of our model. And in Sect. 6.7 we describe some conclusions. 



6.2 The radiative Instability model of Owocki &; Rybicki 

In this section we outline briefly the model for the radiative instability of Owocki & Rybicki J52|, which wc 
are going to use for our linear analysis of radiative instabilities in magnetic winds. A more detailed treatment 
including several refinements can be found in the original paper and in the subsequent series of pap ers |53| , 
p>4] , |63|, p5| . We start our analysis with the radiative acceleration by a single line as described by Eq. 3.19 



Sect. |3.3| . If we assume that the average wind is perturbed by a small, linear fluctuation, we can describe the 
perturbed line forceQ — gh + $gL to first order in the perturbed optical depth r = f + 5t and the line profile 
4>{x) — 4>{x) — 4>' (x)Sv r (r)/vth perturbed by the Doppler shift due to the velocity fluctuation. But we assume 



Variables marked with a refer in this chapter to the unperturbed wind. 
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that the temperature remains constant. 

*flL(r,KL) = -0thin / (<t>(x)5T{x,r) + <f>Xx)??^) e- T ^' r Ux (6.1) 

J-oo V V th J 

Now we can integrate by parts the second term of the integrand using the Sobolev approximation for the 
quantities of the unperturbed wind. 

Sg L (r, kl) « -ffthin / <p{x) (st(x, r) - e^^ dx (6.2) 

We assumed that the line opacity kl is not affected by the perturbation. From Eq. 3.2C we can derive the 
perturbation of the optical depth 

<5r(i,r) = / [ K\ J 8p(r)4>{x) — p(f)<j>'(x) — - — ) df (6-3) 
Jr. V w th / 

V P v A r ) J 

/>oo 

~ Tl / ^—(j)(x) dx (6-5) 

Jx T L 

Here we used the same trick as for 8g^ to get rid of the </>' term. In the WKB approximation, that gradients of 
the mean flow are small, we can use 

Stl _ Sp _ Sv' r (f) _ Sv' r (r) ^ 
tl p v' r (r) v' r (f) 

In our linear analysis we can describe the perturbation by a sinusoidal wave 6v r (r) ~ bv T e %kr . In the limit of a 
strong line (tl 3> 1) we can find a simple expression for Sg^. 



R 



SgL _ ikL s g thin T L 

Sv r Vth 

As in Chap. ^ we can assume that i is close to x and therefore use 



x)e- Ti *< r) / (f>(i)e lkf dxdx. (6.7) 



v' (r) 

x ~ $ - !?rVl(r - r). (6.8) 



This leads to 

$9l _ ikL s g thin T h 
Sv r (r) v t h 



(p(x) 

OO 



)e- T(£ - r) / 4>(x)e- lkLs{i -^ dxdx. (6.9) 



For a strong line (tl ^> 1) Eq. 6.9 can be approximated analytically. We shift the integration from frequencies 
(x) to optical depths (r) using 

dr(x,r) 

— ^- — = -TL<f>{x). (6.10) 

The exp r(x,r) term restricts any significant contribution to the integral to the frequency range where f — f{x,r) 
is close to one. So we define x b by 

f(x b ,r) = l. (6.11) 

We can now express x and x to first order by 

x-x = -^-^- (6.12) 

TL0(X 6 ) 

This leads to 

Sg L ikL s g thin f° _ f f° ( f - f \ 

- — -- = / e / exp [ikL s — — — rfrrfr (6.13) 

0v r [r) w th r L J fL J f \ (p{x b )T h J 

9mM f° exp (_ ( Jk^ e -f d , (6.14) 
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where we used the convention of Owocki & Rybicki [f52f 



ffthin^(a;fe) 

L0 b — 

X h = pn L (j)(x b ). 



(6.16) 
(6.17) 



The radiative acceleration of a hot star wind is not properly described by the acceleration due to a single line. 
Therefore we can not expect that Eq. 6.15 properly describes the perturbation in the radiative acceleration. It 



is rather necessary to sum Sg^ for an ensemble of lines using the same recipe as CAK used for the unperturbed 
wind. The CAK model takes the large number of weak lines into account, as well. Therefore we can not use 



Eq. 5.15, when we integrate over the line distribution function. But for a Doppler line profile 



4>v{x) = —=e~ 



Eq. 5.9 can be integrated to 



ikL s g th[n T L 



5v r (r) 2v t hV^ 
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(6.18) 

(6.19) 
(6.20) 



where 

2 f°° _ 2 
erfc(x) = —= I e y dy 

V n Jx 

is the complementary error function. We find then for the perturbed line acceleration analogous to Eq. 3.33 
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(6.22) 



(6.23) 
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Owocki & Rybicki j52| analyzed this integral numerically. From this analysis they found that Eq. 6.24 can nicely 
be approximated by the analytical expression 



5v r (r) 



= n 



ik 
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(6.25) 



which is just a modified version of the analytical expression for a single strong line (Eq. |6.15 ). The combined 
amplification rate of all lines is given by 



^OR = 
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(6.26) 



(6.27) 



is the mean blue-edge absorption strength. Here c is not the speed of light but an empirical parameter: c « 1.6. 

Owocki & Rybicki Q were criticized by Lucy Q for the neglect of damping due to diffuse radiation. But 
Owocki & Rybicki |^3| showed that this effect reduces the amplification rate only by approximately 50% at one 
stellar radius and 20% at infinity. In this initial analysis we emphasize the effect of the magnetic field. Therefore 
we choose the simple description of Owocki & Rybicki f32| instead of the more exact but also more involved 
description of Owocki & Rybicki |53| or the description of Gayley & Owocki p2| , who analyzed the instability 
in optically thick winds. 
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6.3 The dispersion relation 

To analyze waves in the wind of hot stars we start from the equations of magnetohydrodynamics for a com- 
pressible, nonviscous, perfectly conducting fluid as described in Jackson p7| and add the analytic description of 
Owocki & Rybicki |52] for the influence of the stellar radiation on the plasma. We start with the time dependent 
equations for a nonviscous, perfectly conducting gas: 



9 ~dt + P(V ' V)V + Vp 
dB 

~dt 



B 

47T 



x (V x B) - F rad 



V x (v x B) 



(6.28) 
(6.29) 
(6.30) 



where F ra( j is the force term due to radiation pressure. Now we express the pressure p by v 2 p, where v s is the 
speed of sound, and replace p, v, and B by sums of a equilibrium value and a perturbation. In the comoving 
reference frame (v = 0) this is: 



P = P + M x >*) 

v = Sv(x,t) 

B = B + 5B(x,i) 



Equations 6.28- |6.30| give to first order in small quantities: 

dSp(x, t) 



= 



= p 



= 



at 

<9<Sv(x, t) 

dt 
— <^F ra d 
dSB(x,t) 

at 



+ pV • <5v(x, t) 



B 



+ v;V5pCx, t) H x (V x <5B(x, t)) 

in 



- V x (5v(x,t) x B) 



(6.31) 
(6.32) 
(6.33) 



(6.34) 

(6.35) 
(6.36) 



We assume now, that the radiative force acts only in radial direction and depends only on the radial velocity. 
These equations can be reduced to a single equation for 5v: 







<9 2 £v(x, t) 
dt 2 



w s 2 V(V-(5v(x,i)) + 



va x V x [V x (<$v(x, t) x va)] 

1 dSv r (x,t) 

p dt dv r r 



(6.37) 



with the vectorial Alfven velocity va = B/y/Anp. If we use the result of Owocki & Rybicki |52| for the linear 
perturbation of the radiative force, we can get a dispersion relation from this wave equation for plane waves: 



k-x— iujt 



Equation 6.37 then becomes: 





5v(x, t) = Sve 1 



-lj 2 5v + (u s 2 + u A )(k • <5v)k + (v A • k)[(v A ■ k)5v 



(6.38) 



-(va • <5v)k — (k • <5v)va] 



X OR + 



7j— w((5v ■ e r )e r 



(6.39) 



Equation 6.39 is a vector equation, which is linear in Sv. We can think of it as a generalized eigenvalue 
problem: 

(A(k, v A , v 2 s ) - coB(k r , O R, X nR ) - ^ 2 l)^v = 0, (6.40) 



where A and B are tensors. Wc can find u> and <5v numerically. Then SB and 5p follow from Eqs. 6.34 & 6.36 

P 



dp 
SB 



-(k-5v) 
[(k-(5v)B- (B-k)(5v] 



(6.41) 
(6.42) 
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Although Eq. 6.40 is very involved in the general case, we can find an analytical solution for a simplified situation: 



k(r) = ke r 
v A (r) = ^Ae^ 

The latter approximation is quite accurate far away from the star. In this limit we find 



to 



± 



^ORfc 



+ (Vs + v\)k 2 . 



2(x OR +ifc) V V2(x OR +ifc) 
In the long wavelength limit (k <C X OR ), where the waves are stable, this leads to 



(6.43) 
(6.44) 



(6.45) 



(6.46) 



In the case of a weak magnetic field with v r 3> v s > T>a, this resembles Abbott's Q| result for stable radiative- 
acoustic waves with a fast inward and a slow outward mode. In the case of a strong magnetic field we have 
v r « T>a 3> v s . This leads to higher phase velocities for both modes and reduces the relative difference between 
inward and outward waves. Additionally inward running waves are not advected outward by the average wind 
motion any more, because their phase velocity is higher than the velocity of the unperturbed wind. In the short 
wavelength limit (k 3> Xqr ) we ^ n( ^ 



. ^or , 

LO = 1 ± 



.% + (t£ + «i)*2. 



(6.47) 



These waves propagate, if ^or is less than 2k\J v% + v A , with the same phase velocity inward and outward. The 
amplification rate is also the same for both modes. 



In the case of no magnetic field Eq. p. 40 leads to 
= lu 3 



OR 



ik r 



2 2, 2 

oj — v a k u> 



9-c 



X, 



OR 



jjTVs(kl + k 2 e ). 



(6.48) 



For kg = k^ = this reproduces the result for isothermal waves found by Owocki & Rybicki |5^ |. For oblique 
waves there is a third wave mode. 



6.4 The unperturbed wind models 

To analyze the effect of these waves we construct three wind models for a standard massive star. Our model has 
M = 23M Q , L = l.7x 10 5 L Q , R = 8.5i? Q , T = 60000K, a ca k = 0.56, and fc ca k = 0.28. The first model (model A) 
is the standard analytic CAK wind for a star without magnetic field, rotation, or pressure. In this model we 
reproduce the previous result of dominantly inward running waves found by Owocki et al. In the second 

model (model B) we add a radial magnetic field. Since this magnetic field is parallel to the natural stream lines 
of the wind of a nonrotating star, this field does not change the velocity profile of the unperturbed wind. But it 
changes the microphysics for waves. The last model (model C) is a luminous fast magnetic rotator model. The 
rotating magnetic field provides an additional driving force, which changes the properties of the wind drastically. 
The terminal velocity decreases and the mass loss rate increases. The wind efficiency (Mw 00 )/(L/c) is 3.8, which 
is much higher than for a purely radiatively driven wind in the single scattering limit, where the efficiency can 
not be higher than unity. In spite of the known spin-down problem we choose a strong magnetic field and a high 
rotation rate in order to emphasize the influence of these parameters. For this model we use the model from 
Chap. | (Eq. pHS^). Table |0 gives the magnetic field, rotation rate, and the resulting values for the above 



mentioned unperturbed wind models. In Sect. 6.6 we will discuss how these values and their observation can be 
influenced by waves. In this work we do not discuss the influence of the waves back on the unperturbed wind 
as Koninx |3fj] did. A model with rotation but without magnetic field is not included, because this model has 
the same microphysics for waves as model A. Just M and the velocity dependence on r are different. Model C 
differs from model B in the local conditions by the fact that, due to the rotational twist, the magnetic field and 
the direction of the radiative force are not parallel anymore. Even at the base of the wind is approximately 
-1.67B r . 
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Table 6.1: Results for the unperturbed wind models 



Model 


A 


B 


C 


B r o [G] 





500 


500 










0.92 


Woo [km s _1 ] 


1146 


1146 


784 


M [10- 6 Af o yr- 1 ] 


0.6 


0.6 


17 


R T =2/3 [Rq] 


8.5 


8.5 


11 


(M Voo )/(L/c) 


0.2 


0.2 


3.8 



6.5 Numerical results for waves 

Since our unperturbed wind model is limited to the equatorial plane we limit our discussion for waves to the 
same plane. We have still two free parameters then: The wavelength and the azimuthal angle of k. First we 
will discuss radial waves. This will show all important properties of this wave model. At the end of this section 



we will briefly discuss the influence of the azimuth angle of k. Figures 6.1fe |6.2| show numerical results for 
wind model A&B. The sonic wave modes^j found for model A are identical to the slow magnetosonic modes of 
model B. For model B we find additionally four fast wave modes. These modes, the fast magnetosonic and 
the Alfvenic modes, have the same phase velocity vph ~ va, because both, B and k, are parallel to e r . But 
these wave modes are stable and show no dependence on wavelength. Therefore we concentrate our discussion 



for model A&B on the slow waves: Fig. 3.1 shows the dependence of the waves on the wavelength for r = 2_R, 
plotted relative to the Sobolev length, which is the relevant length scale for radiative wave amplification in the 
model of Owocki & Rybicki |5^] , we use here. R — 8.5i?© is the stellar radius. In the long wavelength limit 
we find stable waves with a high phase velocity inward and a low phase velocity outward. This reproduces 
Abbott's [Q result of radiative-acoustic waves. In the short wavelength limit we find the same amplification 
timescales and phase velocities (except the direction) for both modes. In this limit we would expect to see both 
wave modes in the wind. But the most interesting case is the bridging case A ~ Lg. Here we find the shortest 
amplification timescale for all wavelengths. Inward waves are two orders of magnitude faster in amplification 
than outward waves. They have also a higher phase velocity. This resembles the result of Owocki et al. pt| , 
who did a nonlinear calculation and found that primarily the inward running waves steepen into reverse shocks, 
which are advected outward. The aim of this chapter is to argue that this scenario changes if a magnetic field 



and rotation are involved. Figure 6.2 shows the radial dependence of the wave with A = Ls(r) = vth/v' r . The 
phase velocity for the slow magnetosonic modes is approximately the velocity of sound, so that inward running 
waves are advected outward in the supersonic part of the wind. The amplification of the waves is strongest close 
to the star, where the Sobolev length is short. The phase velocity of the fast modes goes with r _1 for large radii 
since B = B r ~ r~ 2 . This might lead to a small velocity for outward running shocks at large radii. 



Figures p.3p z6.4 show the same plots for model C - with magnetic field and rotation. The crucial point 
is that the magnetic field and the amplifying stellar radiation are not parallel anymore. For large radii they 
are even perpendicular. Therefore the fast magnetosonic wave modes, which are most interesting for us, are 



amplified as well. Figure 6.3.b shows six modes with different phase velocities. Two of them, the Alfvenic 



modes, show no dependence on wavelength. They are unaffected by the radiation field and therefore stable. 



Figure O.a shows the amplification timescales for the magnetosonic waves. The fast magnetosonic waves grow 
approximately one order of magnitude faster than the slow waves. We can therefore expect that these waves will 
dominate. The crucial point is that they are much faster than the unperturbed wind - especially close to the 



star. Inward running fast waves will therefore not be advected away from the star. Figure 6A.b shows that the 
phase velocity for the fast magnetosonic waves remain high for large radii. This velocity is a lower limit for the 
velocity of outward running shocks. Rybicki et al. |63| showed that non-radial perturbations in the stellar wind 
are damped close to the wind. But the magnetic field of our model C is mostly tangential already at the stellar 
surface. For the fast magnetosonic modes |<5i^/<5iv| is 0.5 at the stellar surface and 0.35 at r — 2R. We expect 
therefore, that the effects found by Rybicki et al. will influence but not completely dampen the magnetosonic 
waves. We emphasize that fast magnetosonic waves propagate fastest perpendicular to the magnetic field with 

It is very speculative to draw conclusions for nonlinear waves and shocks from a linear stability analysis. 



2 Ordered by the radial phase velocity we denote in our figures the fast magnetosonic modes with 1&6, the Alfvenic modes with 
2&5, and the slow magnetosonic modes with 3&4. The sonic modes of model A, which are identical to the slow magnetosonic modes 
of model B, are denoted with 3&4 as well. 
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Figure 6.1: Amplification timcscalcs (a) and phase velocities (in the frame of the unperturbed wind) (b) versus 
wavelength for model A&B and r = 2R. In model A only modes 3&4 exist. Modes missing in Fig. a) are stable. 
Mode 3 with A s=s L$ has the shortest amplification timescale and therefore will dominate the wind. These inward 
running waves are advected outward with v r = 811km s _1 and steepen into reverse shocks. 
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Figure 6.2: Amplification timcscalcs (a) and phase velocities (in the frame of the unperturbed wind) (b) versus 
radius for model A&B and A = is- In model A only modes 3&4 exist. Modes missing in Fig. a) are stable. 
The wave amplification is strongest close to the star, where the radiation field and wind acceleration are strong. 
Inward running waves originating there will be advected outward and steepen into reverse shocks. 
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Figure 6.3: Amplification timcscales (a) and phase velocities (in the frame of the unperturbed wind) (b) versus 
wavelength for model C and r = 2R. Modes missing in Fig. a) are stable. In the long wavelength limit all modes 
are stable. In the short wavelength wavelength limit the fast magnetosonic modes (1&6) grow approximately 
one order of magnitude faster as the slow magnetosonic modes (3&4). We expect that fast magnetosonic waves 
in both directions are equally dominant in the wind of model C. Since the magnetosonic waves are much faster 
than the unperturbed wind (v r — 294km s _1 ), the inward running waves will not be advected outward. The 
stable Alfvenic modes (2&5) show no dependence on wavelength. 
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Figure 6.4: Amplification timescales (a) and phase velocities (in the frame of the unperturbed wind) (b) versus 
radius for model C and A = O.lLs- The Alfvenic modes (2&5) missing in Fig. a) are stable. The wave amplifi- 
cation is strongest close to the star, where radiation field and wind acceleration are strong. Fast magnetosonic 
waves originating here will dominate the wind. They have the shortest amplification timescale for all radii and 
a high phase velocity even at large radii. This may lead to fast shocks running outward, which have a strong 
influence on observation. 
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But we showed that a magnetic field has a strong influence on the linear analysis. Therefore we expect a strong 
influence of the magnetic field on nonlinear perturbations as well. Our speculative scenario for waves in the wind 
of hot stars is the following: Waves are predominantly generated close to the star where the stellar radiation 
field is strong and the Sobolev length is short. Waves with a wavelength short compared to the Sobolev length 
are generated predominantly, because they have the shortest amplification timescale. Inward running waves run 
into the star and disappear because their phase velocity is higher than the unperturbed wind velocity. Outward 
running waves, which have, for short wavelengths, the same growth timescale than inward running waves, can 
run over many stellar radii and grow. For the outward running waves Sp and Sv r are in phase. Therefore they 
can steepen into forward shocks with regions of high density at high velocity and influence the observation of 
the terminal velocity and mass loss rate. 

The picture drawn in the last section does not change significantly when non-radial waves are taken into 
account. There are two new effects compared to the radial case. A third wave mode appears in model A as 
predicted by Eq. 6.48. This inward running mode has a long amplification timescale and a relatively low phase 
velocity. Therefore we do not expect a significant contribution of this mode to the situation in the wind of 
model A. For model B we find that Alfven and fast magnetosonic waves have different phase velocities now. 
The fast magnetosonic waves are now amplified, too. But their amplification timescale is much longer than for 
slow magnetosonic waves. 



6.6 Consequences of our model 

In this chapter we do a linear stability analysis and show that the waves in the stellar wind can help to understand 
the observations. In order to calculate quantitative results, which can be compared with observations, it would 
be necessary to analyze the detailed properties of the shocks resulting from the waves found in this chapter. But 
we can derive some estimates from our calculations. 

From Fig. |6.4[ .b we see that the phase velocity of the outward running waves and possibly shocks in the rest 
frame of the star is at least twice the terminal velocity of the unperturbed wind. In the outward running waves 
the oscillations of v r and p are in phase. Therefore a significant amount of matter, but presumably not all matter 
will escape at this or a higher velocity, if the outward running waves steepen into forward shocks. The terminal 
velocity will then be overestimated at least by a factor of two in the observation. For our model C this would be 
^oo.obs ~ 1500km s" 1 . Krolik & Raymond |3l| found that in a nonmagnetic wind shock shells are running much 
faster than the phase velocity of the waves. In an unperturbed magnetic wind model such a high value for Voo 
combined with an reasonable high value for M can only be obtained with a very high magnetic field and fast 
rotation, which leads to a spin-down problem. 

The influence of our model on M is more difficult to estimate. We gain a real factor on 28 in M between 
our model A&B and model C even in the unperturbed wind due to the driving force of the rotating magnetic 
field. Furthermore the observation of M is influenced by the clumping of the wind matter. But to calculate the 
clumping factor <p 2 > / <p> 2 at large radii, where radio observations are made, it would be necessary to do a 
nonlinear calculation including large radii, because the waves steepen into shocks very rapidly due to the short 
amplification timescales. This is beyond the linear model presented here. 



6.7 Conclusions 

In this chapter we analyzed the interaction between a magnetic field and linear waves induced by the radiative 
instability. We found both models complement each other. The magnetic field suppresses the inward running 
waves, which dominate in nonmagnetic winds. This may allow the outward running waves to support the 
unperturbed wind as described by Koninx ]30] | and to form high density shock shells running out at a high 
speed. These shock shells may explain the high terminal velocities measured in winds of massive stars. The 
outward running waves will also lead to an overestimation of M due to wind clumping and the overestimated 
Voo- Wind clumping also occurs without a magnetic field. But in this case the resulting shock shells will run 
inward; and the argument about the overestimated would not apply. The overestimated M and Voa put 
unnecessarily strong restrictions on fast magnetic rotator wind models. We argued that even for a magnetic 
field with B r o = 500G the spin-down time is consistent with the lifetime of the star inferred from the mass loss 
rate considering the uncertainties in the stellar structure. From the observation of nonthermal radio emission 
in many OB and Wolf-Rayet stars we know that these stars have a non-negligible magnetic field. Further direct 
observations are necessary to infer the actual strength of these fields. Previous observations using the Zeeman 
effect were due to the strong Doppler broadening of the lines not sensitive enough to measure reliably the 
magnetic field strength in the winds of O and Wolf-Rayet sta rs pgL It might be possible in the near future to 
measure these magnetic field strengths using the Hanle effect |25|, |26f| . 
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Figure 6.5: Expansion time versus radius. Model C has a slower acceleration close to the star and a lower 



terminal velocity. A fast magnetosonic wave with an amplification timescale of about 2500s (cf. Fig 3. 4. a) grows 
by a factor of » e 40 , while the wind expands from r = 1R to r = 1.5R. The electron scattering opacity is 2 jz at 
r = 1.3i?. Therefore we can expect shocks already at this radius, where line observation starts. The objections 
of Lucy pTj do not change this situation qualitatively. 



We showed that a luminous fast magnetic rotator model plus wind perturbations by waves or shocks can help 
to explain the observed high values for M and i>oo without being ruled out by the spin-down problem. Further 
observation of the magnetic field and further theoretical work on the evolution of stellar rotation are necessary 
to evaluate to role of magnetic fields in winds of massive hot stars. 



Chapter 7 

The two dimensional fluxsheet model 



7.1 Introduction 

In Chap. | we discussed an equatorial wind driven by radiation and magnetic fields. In this chapter we extend 
our discussion to the non-equatorial region. This is necessary to get better quantitative results for mass and 



angular momentum loss (cf. Sect. 4.5). It has also been argued that the observed high mass loss rates and 
terminal velocities are not created in the same regions of the wind |6l], ^2). Furthermore it is possible to apply 
a non-equatorial model to winds from accretion disks. This can be important for some types of stars and for 
Active Galactic nuclei (AGN). It is obvious that disks rotate close to the critical (Keplerian) rotation rate. And 
it is widely assumed that the magnetic field is a major mechanism to extract angular momentum from the disk 
and to drive winds and jets jll} . But the role of radiation from the central object or from a hot disk is often 
underestimated. 

In this chapter we derive a model for a non-equatorial wind driven by radiation and magnetic fields. Our 
model is again based on the fluxsheet concept and therefore essentially one dimensional. But we generalize the 
model of Chaps. ||&|| by dropping the restriction that the fluxsheet has to lie in the equatorial plane. In App. |A| 
we will show the equations of Chaps. follow from the generalized equations developed here. At the end of 
the chapter we will show results from first numerical exploration of the model. 

In the first part we follow closely the derivation of Lovelace et al. |35| to reduce the equations for a rotational 
symmetric, rotating, stationary, and perfectly conducting MHD flow to the fluxsheet concept. Since we are 
primarily analyzing winds from stars and not from disks we shifted the derivation from cylindrical to spherical 
coordinates. The goal of Lovelace et al. is to derive an ordinary differential equation for the flow in the fluxsheet 
from Bcrnoullis equation of energy conservation. This is reasonable for a pure MHD model. But we can not go 
this way because, opposite to Lovelace et al., we want to include radiation as a major wind driving mechanism. 
We saw in Chap. | that we have no simple theory for the energy exchange between radiation an matter. The 
Thompson and the CAK theory, both give us only a source term for momentum. Therefore we leave the path 
of Lovelaces model and derive a final ordinary differential equation from Eulers equation of momentum balance. 
This equation is a vector equation. We take only the component parallel to the a priori given fluxsheet into 
account. The component perpendicular to the fluxsheet can be used to improve the shape of the fluxsheet. 

In Sect. [7.2| we follow to derivation of Lovelace et al. for the basic equations in the fluxsheet description. 



In Sect. [7.3| we introduce our fluxsheet coordinate system with its a priori given shape. In Sect. 7_A we use the 
results of the two previous sections to derive our final wind equation as a single ordinary differential equation. In 
Sect. (7^ we derive the transverse component of the Euler equation, which was neglected in the previous section. 
It allows an improvement of the chosen fluxsheet shape, and finally we will give some preliminary numerical 
results in Sect. |7.6| to show the capabilities of this model. 



7.2 The basic model equations 

As in Ch ap. fl| we start from the equations of ideal magneto-hydrodynamics (MHD) in a quasi stationary situation 



(Eq. 4.1-4.5). Again we assume that the wind is rotational symmetric and neglect therefore all dependencies on 
the azimuthal angle cf>. But now we do not set vg and Bg equal to zero as we did in Chap. ^. Rather we use the 
perfect conductivity equation 

E+-vxB = (7.1) 

c 
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and the axisymmetry condition E$ — to find 



v r e r + vge e 
K(r,0)B p , 



(7.2) 
(7.3) 



where the subscript "p" denotes the poloidal component of the vector. From = V ■ (pv) = V ■ (pv p ) = 
V • (pn(r, 6*)B P ), we know that 

B P 'VHM))=0, (7.4) 
because V • B p equals V • B = 0. The poloidal magnetic field can be expressed in terms of the flux function 

f = r sin 9 A^r, 9), (7.5) 

where Aa, is the azimuthal component of the vector potential. ^ is constant on every fluxsheet, i.e. 



From Eq. 7.4 we get 



B p • V* = 0. 



(7.6) 



(7.7) 



where F(^f) is a fu nction whose only relevant feature is that F is constant along a fluxsheet. Using \V^f\/(r sin 
= B p and Eq. we get 

v x B = — !— -(t^, - kB )V*. 



r sm ( 



From Eqs. [T|&]7j| we know tnat 



r sm 



j(«0-/tB^)=w(*), 



where u;(\&) is another arbitrary function of 'J'. From Eq. 7.8 we get the electric field 



E = -iw(*)V* 

c 



We can now express the velocity field by 



v - ^ Br 



Anp 
Airp 



inp 



B<f, + rsin0u(^) 



B + <*;(*) x r, 



(7.8) 
(7.9) 

(7.10) 

(7.11) 
(7.12) 



where u> equals Lu( i S)e z . And e z is the unit vector parallel to the rotation axis. At the base of the wind, where 
p is large, va, equals uirsinO. Therefore Q is the angular velocity of the rotating fluxsheet. Using the toroidal 
component of Eq. 4.2 we find 

B p -V(rsin0(B -Ft; )) = O, (7.13) 



which implies that 



H(V)=rsin9(B <f> -Fv 4> ) 



(7.14) 



is another arbitrary function of This corresponds to the conservation of angular momentum for matter moving 
on a given flux surface, = const. Now we can express va, and Ba, by 



r sin 9Ba, = 
r sin Ova, = 



H + r 2 sin 2 6 Flo 



1 - Ml 



r sin Olj 



FH\ 
Anp ) 



where we have used the poloidal Alfven Mach number 



/47rpv| 



F 



(7.15) 
(7.16) 

(7.17) 
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On the Alfven surface (r = r Aci — #Ac, M p — 1) the denominator of Eqs. 7.15 & [7.1q vanishes. Since B$ and 
v<p must be finite on the Alfven surface, we get the conditions 



H 



r\ c sin 2 9 Ac u) 



-r Ac sin 2 9 Ac Flu 
FH 



From l'Hospital's rule we find 

(rsin^)A = -H 



47rpAc 
d[ln(r 2 sin 2 6»)] A /ds 



d{hxp) A /ds 



(r sin8v^)A = (rsm9)\uj 
= (rsm9) A uj 



1 + 



d [ln(r 2 sin 2 9)] fds 
d(ln p)x/ds 



_ (r sm9Bj,) A 
H 



where s is the arc length along the poloidal projection of the field line, i.e. d^f/ds = 0. 



(7.19) 
(7.20) 

(7.21) 
(7.22) 
(7.23) 



7.3 The fluxsheet coordinate system 

Up to this point our model is exact besides the approximation of rotational symmetry, stationarity, no viscosity, 
and perfect conductivity. But therefore we have still a system of two partial differential equations with two 
dependent (v r , v$) and two independent (r, 9) variables. Such a system was solved by Sakurai for a stellar wind 
without radiation pif . 

In order to reduce the problem from two to one dimension we split the three dimensional space into a one- 



parameter family of two dimensional fluxsheets. Figure 7.1 shows how these fluxshects fit into each other like the 



leaves of a blossom. Every fluid element will flow along a single fluxsheet. The magnetic field lines are confined 
in the individual fluxsheets as well, because we assume flux freezing and a quasi stationary model. Due to the 
rotational symmetry of the whole wind, the fluxsheets are rotational symmetric as well. Fluxsheets exist even 
in a model which solves the full MHD equations. In this case we would find simultaneously the shape of the 
fluxsheets and the motion of the matter within the fluxsheets. We already reduced the problem from three to 
two spatial dimensions by assuming rotational symmetry. Now we reduce it again by one spatial dimension by 
using an a priori given shape for the fluxsheets and only asking for the motion within the fluxsheets. Falcke JlS| 
integrated the equations of an radiation driven nonmagnetic wind simultaneously for the wind velocity and the 
shape of the fluxsheet. This is much more difficult for a magnetic wind, since our solution has to pass through 
three instead of one critical point. And this will only happen once we found the two eigenvalues of the equations. 
The wind solution has to pass through three critical point if it should extend from the subsonic base of the wind 
to arbitrary large radii. 

We require that the shape of the fluxsheet is given by r(u, rj) and 9(u,r]), where rj is constant along the 
fluxsheet and u is an arbitrary parameter along the fluxsheet. S ince the parameter r\ is constant along a 
fluxsheet, it depends only on the flux function \I/ introduced in Eq. (7^ and vice versa. The reason for using rj 
instead of the already introduced parameter \1/ is that we can choose it arbitrary in the same way as u. It just 
distinguishes different fluxshects, while \E' has a deeper physical meaning. For our purpose it is not necessary to 
know \& explicitly. Therefore it is numerically simpler to introduce rj. rj and u are a curved coordinate system for 
the meridional (r, 9) plane which is adapted to the shape of the fluxsheets. It is therefore more appropriate for 
solving the wind equations than the (r, 9) system. We do not give explicit formulas for r(u, rj) and 9(u, rj) here, 
so that the shape of the fluxsheets remains yet unspecified. This has the advantage that the wind equations, 
we derive now, remain independent of the explicit shape of the fluxsheet. Later, when we start our numerical 
computations, we have to specify the shape of the fluxsheets through explicit formulas for r(u,rj) and 9(u,rj). 
But then we can easily compute and compare models for different fluxsheet shapes by using different explicit 
formulas for r(u,rj) and 9(u,rj). 

Now we can define a unit vector e p tangential to the fluxsheet in poloidal direction given by 

sv sr9 
e p (u,rj) = , =e. r + —= , e e (7.24) 

V r 2 + r 2 9 2 V » ;2 + r2 8 2 

= e pr (u,rj)e r + e p g(u,rj)e g , (7.25) 

where the dot ' refers to the partial derivative with respect to u. The constant s — \f\/r — ±1 fixes for later 
convenience the orientation of e p so that e p , per definition, always points outward. We can now express v and 



50 



CHAPTER 7. THE TWO DIMENSIONAL FLUXSHEET MODEL 




Figure 7.1: shows the stellar surface and three members of the infinite set of fluxshects. Every fluxsheet has a 
different rotation rate u>(rj). Whether u> is higher at the poles or at the equator is not clear yet. 




Figure 7.2: A non-equatorial fluxsheet in the (r, 8) plane. In this plane the thickness of the sheet is dh{u). The 
rotational symmetry contributes a factor of 27r?-sin0 to the total cross section of the fluxsheet. 



B as functions of the fluxsheet coordinates u and r\: 



B = B p (u,r])e p + B^u.r^e^ 
4irp(u, 77) 



-(v — u(rf) x r) 



(7.26) 
(7.27) 

(7.28) 



In order to describe p and B p in the new coordinates we need an expression for the cross section of a fluxsheet. 
Then we can use the conservation of matter and magnetic flux as we did in Eqs. 4.S&4.9. In Fig. 7.2 we show 
the wind between two infinitesimal separated fluxsheets. The cross section dA between the two fluxsheets is 



dA = 27irsin0|<#i|. 



(7.29) 



The separation between the two fluxsheets dh varies along the fluxshects. But we can relate it to drj, which is 
constant along the fluxsheets: 

drj = || dh(u) (7.30) 
Since 77 is constant along every fluxsheet and dh is measured perpendicular to the fluxsheets, we can express 



7.4. THE ONE DIMENSIONAL EULER EQUATION 



51 



drj/dh by 

^L = \Vr)(r,0)\. (7.31) 
Finally we find for the cross section between the fluxsheets 

dA = 27rr sin 6 (7.32) 
|Vt?| 

Now we use B p dA = const, pv p dA = const, and drj = const between the two fluxsheets to express 

r Ac sin6»Ac « P Ac |V?7| 



rsin9 v p \Vr)\ Ac 



PAc(v) (7-33) 



_ RsinOp \Vr,\ 

Bp r S m9 |V77lo P °- ( } 

The poloidal magnetic field strength at the base of the wind B p q is an input parameter of our models. The matter 
density at the Alfvenic critical point p Ac will be derived from model input parameters later. Alternatively we 
can relate the matter density to the mass loss rate in the fluxsheet 

8M 2nrsm9 

~W = pVp ^T- (7 - 35) 

If we want to calculate the real mass loss rate for our model star we must solve the wind equation for every 
fluxsheets from the equator to the poles and then integrate the values found for dM/drj over all flux sheets: 

^(equator) gj^ 

M rcal = 2 / — dri (7.36) 

•Aj(pole) ar \ 

In order to get a better intuition for the mass loss at different latitudes it is helpful to assume for a moment 
that the density and the wind velocity are constant at the base of the wind. Under this approximation we can 
integrate Eq. 7.3£ neglecting the oblateness of the rotating star and find 

M approx (r?) = 4nR 2 Po (r])v P Q(r 1 ). (7.37) 

This quantity is analogous to mass loss rate used in Chap. [|as one of the eigenvalues we have to find numerically. 



7.4 The one dimensional Euler equation 



The standard way in disk wind physics is to derive v p from Bernoullis equation, which essentially describes the 
conservation of energy in the fluxsheet. We will use Eulers momentum equation (Eq. 4.2) in order to include 
the momentum terms for radiation pressure. Since we have already expressed in terms of u, r], and v p we 
can neglect the 4> component of Eq. 4.2. Theoretically it would be sufficient to integrate either the r or the 9 
component of Eq. 4.2. But this would lead to large numerical errors when the flow direction e p has a significant 



component perpendicular to the chosen component of Eq 4.2. Therefore we multiply Eq. 4.2 with e p and obtain 
a single ordinary differential equation in v p . This equation allows us to integrate v p along the fluxsheet. To 
obtain a solution for the complete wind we have to integrate v p along several fluxsheets. 

Although we neglect for a while the force balance perpendicular to the fluxsheet, we can not neglect that v p 
depends on u and r\. We solve this problem by separating the variables: 



V p (u,T)) = V(r))v(u), 



(7.38) 



where we use an a priori given function V{rj). This is the analogy of the self similar solutions of Blandford & 
Payne, who a priori assumed certain powcrlaw dependencies of p and B on the cylindrical radius. Now we have 
a single differential equation with v(u) as the only dependent variable. 

But this equation still contains the differential operators d r and dg. We transform these into our fluxsheet 
coordinates by 



(7.39) 
(7.40) 
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dr 


= u 7T " 
\ou 
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where we take into account, that p, v^, and B depend on u, 77, and y. For brevity we will use the following 
convention for partial derivatives: 



dr 1 

d ± = f 
du 



dj_ 

dx 



= f 



df 
80 
d£ 
ch] 
df 
dy 



= f 
= f 
= f- 



(7.41) 



Here we used the already the dimensionless quantities x 
multiplied by e p : 



= e p - 



( V . V) V + -5- + V$ - g rad 
P 



r/r* and y = v/v*. We can now express Eq. 4.2 
(V x B) 



— V 



(7.42) 



Doing all these substitutions and simplifying the result to a readable form includes a lot of complicated 
algebra. We used therefore the software package Mat hem atica [fri[ to do most of the steps. With a small 
Mathematica program it is possible to show that Eq. 7.42| is equivalent to 



with 



B i(u,y,y) 



Ai(u, v)v + Bi(u, v,v) = 
4irp _ 



(ru' + Ou) 

(fu' + Ou) 

(ft]' + f]9) 
4irp _ 



vV 1 



V 



F 2 



■V<f, 



dv v 



(7.43) 
(7.44) 



, 2 P 



- + d u v 2 
P 



F 2 



VVv 2 



v<t> [Vj, — cur sin 9 + v<j, 



Anp 






(- 

VP 


"p v\ 




,P v) 


+ < 




V 


(p - 


V 



47rp_ . _ 2 

-p^v^ {rsmOuj -v^,) +v rj> 



tan 6* 

f(V1>-g rad (u,v,v)) (7.45) 
— r sin 9ui. (7.46) 

We can make this equation dimensionless by multiplying it with l/{x 2 VQ) — r 2 /{GM) and introducing the 
length and velocity scales we used already in Chap. 0. Here we need additionally the special rules 



V<f, 



vV 



v*yY 
r* _ 



where Y(rj) is just a dimensionless version of V(rj). We get then 

A 2 (u,y)y + B 2 (u,y,y) = 



B 2 (u,y,y) 



(xu + Ou) 
(xu + Ou) 
(xfj + Of}) 



y m 2 v y 



(7.47) 
(7.48) 

(7.49) 
(7.50) 



2 P , o 2 1 V<t> I ' 

y s -+d u y s + [ -yt + y,, 



d v y 2 s 



jj^ixsmOu-y^+yl 



x tan( 



x ( -5 - a lad (u,y,y) 



(7.51) 
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If we use the line-driving mechanism of Castor, Abbott, and Klein |l6j plus electron scattering as description 
for the radiative force, a ra d is given by 



a r ad(w,2/, y) 



GM 9rad 

~ 2 



a; 



KThfthP?"* 



<9;F 



r , C(u,y) 

2 



.r 



dy x 



dx 



= = {ye pr + ye pr ) Yu + (Ye pr + Ye pr ) yr). 



(7.52) 
(7.53) 
(7.54) 
(7.55) 



The derivation of the CAK force in Chap. || shows that the CAK force depends on the gradient of the velocity 
in the direction of the driving radiation. We want to analyze the wind from a star. Therefore the photons 
propagate in radial direction. If this model would be applied to systems where a disk dominates the radiation 
field, the equation have to be modified at this point. Outside the equatorial fluxsheet the poloidal flow v p and 
the poloidal gradient d u are in general not in radial direction. This leads to the complicated relation between 
the two spatial derivatives of the wind velocity, \y x \ and y. Using \y x \ instead of y is the simplest approximation 
for the fact that in a curved fluxsheet geometry the radiation passes through diffe rent fl uxshccts. In ord er to 
compare the results from this chapter with the results from Chap. || we rewrite Eq. 7.49| analogous to Eq. 3.3E . 
This allows us to use the same arguments about critical points and reduces the necessary modifications in our 
software. Since the CAK force term is now proportional to |y x | acak instead of y a ^ k we use \y x \ as the local 
unknown quantity. This leads to 



A(u, y) \y x \ + B(u, y) = C(u, y) \y x \ 



(7.56) 



Before we can use this equation we have to transfer A 2 and Bi into A and B. Once we know \y x \, we need a 
formula which gives us y as function of \y x \. This requires to distinguish the cases y x > and y x < 0. Since 
we are only interested in expanding winds in reasonable fluxshccts we can savely assume that the first case is 
always true. Then we find 



A(u,y) 
B{u,y) 



Ao 



ye pr Yu + (Ye pr + Ye pr ) yfj 



B 3 -A 2 - 
y x - ye pr Yu - (Ye pr + Ye pr ) yfj 



e pr Yu 



(7.57) 
(7.58) 
(7.59) 



where Bz{u,y) = B2 + Ci/" cak is just B2 without the line acceleration term, which appears separately on the 
right hand side of Eq. 7.56 . The deaccelerating solution (y x < 0) differs from the result above only by a few 
signs. The auxiliary quantities of Eq. 7.56| are given by 



P 

M = - 
P 

P 

7 = - 
P 

F 

F 



M 2 

P 



V4> 



|Vry| 

dy\Vr]\ 
|Vr/| 

1 PAc 

2 PAc 
PAc 

P 



tan 8 



tan (9 



x 
x 
x 
x 



Y 
Y 



PAc 
PAc 



xsinf? y |Vry| 



Ac 



a; Ac sin 9a c yAc |V»?| 



2/rot ^ P Dmt,Ac 

y rot (l-M p 2) 

2/rot 



- -11 2 
il rot, Ac 



y rot (l-M p 2)2 



V 



2/r 2 ot(l-Mp 2 ) 



"2/rot - 



9M p 2 



MM 



(7.60) 
(7.61) 
(7.62) 
(7.63) 
(7.64) 
(7.65) 
(7.66) 
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V<t> 

y<t> 

C(u,y) 

c 



dMl y 

. C fxysmOY"^ 



dMl ' p 



„ 2?/ r ot,Acyrot.AcMp 

7^p : 



y rot (l-M p 2) 



V* | V77I 



Ac 



KTh^XAc sm V Ac y A cPAc 



(7.67) 
(7.68) 
(7.69) 
(7.70) 



7.5 The transverse force component 

The model developed in this chapter so far allows us to compute the wind in a fluxsheet of predefined shape. 
This is not very satisfying, because if we choose by chance a fluxsheet whose shape differs too much from the real 
shape we will get wind solutions without physical meaning. On the other hand the numerical effort for finding 
the shape of the fluxsheet and the wind solution simultaneously is much greater. In this case we have essentially 
a 2-dimensional problem. The compromise between these two extreme alternatives is to assume a priori a shape 
for the fluxsheet, to calculate the wind solution and finally to use the transverse component of the Euler equation 
to check whether we have chosen a reasonable shape of the fluxsheet. These two steps could even be iterated 
to get a consistent overall solution. Alternatively we could use the tangential and the transverse component of 
the Euler equation to find the shape of the fluxsheet and the wind velocity inside the fluxsheet simultaneously, 
while we integrate from the stellar surface outward. The difficulty to pass through all three critical points on 
the way would make a completely new computer code necessary. And it would not be clear from the beginning, 
whether this concept would work in the numerical praxis. For our initial calculations we choose a less ambitious 
approach by specifying shape of the fluxsheet a priori. This allows us to use our code for the equatorial wind 
with few modifications. Afterwards we use the transverse component of the Euler equation to check, whether 
the chosen shape for the fluxsheet was reasonable. This allows us to explore the raw shape of the real fluxsheet 
and leads to a much higher confidence in our solutions. For this approach we need a dimensionless form of the 
transverse component of the Euler equation. 



We start our derivation from the Euler equation (Eq. 4J2), whose poloidal component (Eq. 7.42 ) we used to 
derive our wind equation. If we choose the right fluxsheet, the transverse component of Eq. 4^ 



= Ag t = e t 



VP 

v ■ V) v + + V$ - g rad 

P 



e t xJL]. ( VxB, 



with 



e t 



- pr 



eo 



(7.71) 



(7.72) 



should be fulfilled, too. But in general we do not know the exact shape of the real fluxsheet. Therefore Ag t will 
not be exactly zero for our solutions. In order to get an estimate how seriously Eq. 7.71 is violated we can compare 
Ag t with the poloidal or transverse component of (v • V)v. Therefore we use now the same Mathematica code, 
we used for Eq. 7.43, to derive dimensionless versions of Ag t and (v • V)v. We find for the dimensionless version 
of Ag t 



Ao t = -3 I A -9t - e P • ((v • V)v)] 

= A 4 (u,y)y + B 4 (u,y,y) + Ab t 



(7.73) 
(7.74) 



with 
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For the tangential and transverse components of (v • V)v we find 
A6 n 



Je p .((v.V)v) 



(6p r 



x tan 



yyY 2 



° P eV 



y 2 YY 



Ab t 



^e t .((v-V)v) 



x tan o x / 



(e p ee pr . - e pr e p e) I e pr .7? + 



) e pr w + 



(7.76) 
(7.77) 

(7.78) 
(7.79) 

(7.80) 



Now we can use |Aat/A6 t | as relative measure how much the transverse component of the Euler equation is 
violated. But in order to estimate how relevant this violation is for the overall solution we should compare 
the transverse component of the Euler equation with the tangential component. This is given by the quantity 
|Aat/A&p|. 



7.6 First numerical models 

In this section we show some initial calculations, which try to develop further some ideas of Poe et al. J6^] 
for winds from Wolf-Rayet stars. They developed a wind model for Wolf-Rayet stars based on the assumption 
that a fast and thin polar wind, driven only by the CAK line force, is responsible for the high observed terminal 
velocities. In the equatorial plane they assumed a thick and slow wind, additionally driven by a rotating magnetic 
field. Here they used a model similar to ours. This wind should be responsible for the high radio flux, which is 
otherwise explained as large spherical mass loss. Although this model is theoretically appealing it was ruled out 
as general concept for Wolf-Rayet stars by observations which showed that these stars in general do not have a 



strong density gradient between the polar and equatorial regions 1 65 , |66| . In this section we want to argue that 
the model of Poe et al. might be modified to fit this observational fact. The basic idea is that the fluxsheets are 
bent towards the poles and not towards the equator. The latter is assumed in the wind compressed zone and 
disk models [Q. These models assume a nonmagnetic but rotating wind. In this case the centrifugal force will 
cause the deviation from radial outflow. But if we introduce a significantly strong magnetic field, the Lorentz 
force will dominate the centrifugal force. From the physics of jets we know that a rotating magnetic field can 
very efficiently collimate polar outflows. Additionally Owocki et al. have shown that the radiation force has a 
strong non-radial component |57j . This will even enhance the effect of the Lorentz force described in this section. 
For simplicity we do not include this radiation effect in our calculation at the current stage. In any case the 
equatorial outflow will automatically be dispersed. This will have several consequences: (1) The wind density 
in the polar region will be enhanced due to the geometrical compression. (2) The terminal wind velocity in the 
polar region will be reduced due to the fact that the polar fluxtube is collimated in the supersonic region. The 
reduced terminal wind velocity increases the polar wind density as well. (3) Exact the opposite happens in the 
equatorial region. There the terminal wind velocity is increased and the density is decreased. This reduces or 
maybe even eliminates the discrepancies between the polar and the equatorial wind predicted by the model of 
Poe et al. and criticized by the observations of Schulte-Ladbeck et al. [[35j [36) . 

For our numerical models we use here the "Wolf-Rayet test model" of Poe et alii. This model star has a 
mass of M = 13M Q and a luminosity of L = 3 x 10 5 L Q . For the equatorial radius we use a prelimnary value 
of R = 5i?0 for the base of the wind. For strong winds (high rotation rate and strong magnetic field) this will 
lead approximately to a radius of 8Rq for r cs = 1 as used by Poe et alii. For the CAK parameter we choose 
«cak = 0.61 and fc ca k = 0.18. The electron scattering opacity for an Wolf-Rayet star «Th = 0.2cm 2 /g due to the 
missing hydrogen in the wind. For the wind temperature we choose a value of T = 20000K. Poe et al. quote a 
mass loss rate of M bs — 3 x 10 -5 M /yr and a terminal velocity of Woo.obs = 2900 km/s as observational results 
for the wind. These values should be directly or indirectly be explained. 
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Table 7.1: Equatorial wind models for the Wolf-Rayet test star without any extra wind compression or dilution 
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Our wind equations contain the necessary terms to include differential rotation. And from our sun we 
know that stars can rotate differential. It might even be that differential rotation is the normal case for stars. 
Nevertheless we limit our discussion here to the case where the rotation rate of the star does not depend on 
latitude (uj = 0) to reduce the number of free parameters in this initial analysis. 

We can not reproduce the numerical results of Poe et al. exactly as we did for the O star models of Friend 
& MacGregor [^o) in Chap. |5[ because for simplicity we do not include the correction for the finite disk of the 
star Jll| Eq. 50] in our equation at the current stage. Additionally Poe et al. do not print their wind equations 
explicitly. Including the finite disk correction will be neccesary later when this model is developed further to 
explain Wolf-Rayet star winds quantitatively. 

In a first step we calculated a set of equatorial wind models to get a first feeling for the dependence of M and 

on the wind parameter a mt and B p0 . These calculations confirmed the results of Poe et al. that in a simple 
equatorial model we need a high rotation rate to explain the high mass loss rates and a very strong magnetic 
field to explain the high terminal velocities. We agree with Poe et al. that such a model would be ruled out by 
the spin-down problem. 

The next step is to check that the magnetic force really collimates the fluxsheets towards the poles. To 
check this we calculated wind models for a set of non-equatorial fluxsheets. For these first calculations we used 



fluxsheets which were not bent (Fig. 7.4). This fluxsheets configuration would be expected for a nonrotating 



star. This fluxsheet configuration is mathematically given by 



x(u,ti) = Vl±ji ( 7. 8 i) 

u 

9{u,rf) = arccot?7. (7.82) 



It should be mentioned that in this special case the wind equation (Eq. 7.56) becomes independent of Y, p, Co 
and F, because xfj + Orj is always zero. But this is not true for the transverse force component. We calculated 11 
models with rj between and 1 covering an lateral angle of 8 between 90 and 45 degrees. We used a magnetic field 



of 1000 G and a rotation rate of a ro t = 0.95. The latter is more extreme than the values given in Tab. 7.1 and 



Fig. [7_J. But the idea of this section is to propose a model with a strong equatorial wind due to rapid rotation. 
The magnetic filed should be moderate, so that we do not have a spin-down problem. But it should be strong 
enough to spread the mainly equatorial mass outflow of the whole hemisphere, so that this model is not ruled out 
by the lack of observed polarisation. As reference radius for the terminal velocity we chose 100i?. For comparison 
with observation this reference radius should not be chosen too large, because the P-Cygni profiles are not too 
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Figure 7.3: shows M versus Vqq for the models from Tab. 7.1. The solid lines connect models with a constant 
value of a Iot . The dotted lines connect models with a constant magnetic field strength. Additionally the observed 
(+) values and the values for the nonrotating, nonmagnetic case (CAK) are shown. Rotation basically increases 
the mass loss rate. While the magnetic field basically increases the terminal velocity rate. 
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Figure 7.4: shows the meridional plane (r, 8 or u, rj). The fluxsheets start from the stellar surface and extend 
straight outwards. This fluxsheet geometry was used for our first calculations which should prove that the 
fluxsheets are actually collimated towards the poles. The star is oblate due to its rotation (a rot = 0.95). 
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Table 7.2: Wind solutions with three critical points for the fluxsheets shown in Fig. 7.4 
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Figure 7.5: shows Aa t as defined in Eq. 7.73 for the r\ — 1 solution form Tab. 7.2 



far away from the star. The next important step is to check whether the fluxsheets should now be bent towards 
the polar axis or to war ds the equatorial plane. For the r\ = 1 (9 = 45°) solution we evaluated Aat, A6 t , and A6 P 
as defined in Sect. 7.5. The function Y(t]) is fitted to the wind velocities in the different fluxsheets at r — 2R. 
The result is plotted in Figs. 7.5-7.7. The singularity we see in Fig. 7.5 is due to the Alfvenic point, where 
we yet do not use all side conditions to keep the transverse force finite. If we mentally subtract the singularity 
from the plot, we see that Aa t is negative everywhere, and that we get the strongest transverse force component 
close to the stellar surface. In Fig. |7.6| we see the inertial contributions to the transverse Euler equations. Since 
the fluxsheet is not bent we have only the transverse component of the centrifugal force here, which of course 
is directed towards the equatorial plane. The centrifugal contribution has the opposite sign than the magnetic 
contribution from Fig. 7.5. And it is close to the stellar surface more than an order of magnitude smaller. The 
same is true for the poloidal component of the inertial force. We can therefore conclude from here that the 
fluxsheets should be strongly bent towards the polar axis. We find the same result for the rj = 0.1 fluxsheet. 
Although the force towards the polar axis is not so strong close to the equatorial plane. The next step is now to 
calculate wind solutions for fluxsheets which are bent t owa rds the polar axis. As an first approximation for the 
fluxsheets we choose the fluxsheet shape shown in Fig. 7J3. We get again a one parameter family of fluxsheets 
with rj as parameter. The additional parameters A and 5 are constant for all fluxsheets in a stellar wind model. 
The motivation for this nonobvious choice was to produce fluxsheets which are bent towards the polar axis close 
to the star and rather straight far away from the star, where the magnetic field is weak. Our fluxsheets start 
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Figure 7.8: shows geometry of the curved fluxsheets used in this chapter. 

with an inclination of 6 = arccotry for x — > (u 3> S/2) and have an inclination of 77(1 + \-k/2) at x — > 00 
(u <C S/2). So A controls how strong the fluxsheets are bent. While 5 controls how fast they are bent. For A = 
these fluxsheets reduce to the straight fluxsheets used previously in this section. And for A = rj = we obtain 
the equatorial fluxsheet used everywhere in this thesis. The equatorial fluxsheet should always be straight due 
to the symmetry between northern and southern hemisphere. 

We started our analysis of bent fluxsheets trying to find solutions with three critical points. But we found 
that the solutions extend only beyond the outer critical point if the fluxsheets are nearly straight. A numerical 
analysis of the solution around the outer critical points showed that the third condition for the critical point 
(Eq. p3.8[ ) is not fulfilled with the necessary numerical precision. This might have two reasons, (f ) The numerical 
roundoff errors in the formulas for the critical point are too large. Therefore the root finding algorithm of 
our code can not find the position of the outer critical point with the required precision. This would require 
an improvement of our numerical code in order to reduce roundoff errors. (2) There might be no point on 
our fluxsheet where all three conditions for the outer critical point are fulfilled. This would require a deeper 
mathematical analysis of our wind equation, ft might be that we can derive from this problem further conditions 
for the shape of fluxsheets which have a complete wind solution. For this point see also Appendix [b|. Due to 
this problem we calculated few initial wind solutions without the outer critical point analogous to our models in 
Chap. [|. But these models still have fairly straight fluxsheets. The results are still too few and inconsistent to 
be published here. But we showed here that the non-equatorial fluxsheets are indeed bended towards the poles. 
Future numerical computations will complete the picture we sketched here and show how strong the effect of the 
bent fluxsheets really is. 



Chapter 8 

Conclusions 



In this thesis we analyzed new aspects in the theory of magnetic winds from massive stars. This work was 
motivated by the fact that we have many indirect hints for the existence of a significant magnetic field in the 
wind of massive stars although the magnetic field is not accessible for direct observation yet J2q| . The 
development of the theory is still at a very early stage. Therefore we can not yet produce models, which 
allow to describe and understand in detail the complex structure of the winds from massive stars we know 
from observation. But in this thesis we present two new and important components in the theory of magnetic 
winds from massive stars. It is important to recognize that both components are not optional. Both effects are 
automatically present without further conditions, if the star rotates with a significant magnetic field. 

In Chap. ^ we analyzed, in the linear limit, the evolution of wind instabilities in magnetic winds. The 
foundation for this analysis is the well known instability of the line driving mechanism. This instability has been 
intensively studied |5^, |53|, ^4], g5[ Q . And it seems that this instability can explain some details in the spectra 
of massive stars. But in spite of initial hopes |56| no way was found how these waves can help to solve the 
basic wind momentum problem in many of these stars, especially Wolf-Rayet stars. We showed in our analysis 
that the magnetic field is the missing key. It has a strong influence of the behavior of the waves resulting 
from this instability. In the nonmagnetic case the waves are dominantly running inward, leading to reverse 
shocks. If a significant magnetic field is present, the waves are running dominantly outward, leading to forward 
shocks. These forward shocks will contain layers of material which have a much higher velocity than the terminal 
velocity found in stationary and smooth wind models. This is an important contribution to the explanation of 
the observations. Stationary wind models need a very high magnetic field to produce a high terminal velocity 
Such a high magnetic field would cause a spin-down problem, which is easily avoided in our wind models with 
waves. Additionally our waves produce a higher observed mass loss rate by feigning a higher terminal velocity 
and a higher wind density through clumping. This helps to overcome the wind problem in Wolf-Rayet stars as 
well. 

Our second model, the bending of the non-equatorial fluxsheets towards the poles, is a direct and unavoidable 
consequence of significant rotation and magnetic field strength. This effect is a ideal extension to our first model, 
because its results lead into the same direction. We showed that the non-equatorial fluxsheets are bent towards 
the poles. This has several consequences. Due to the extension of the flux cross section in the equatorial region 
we will get an even higher terminal velocity, because the wind is supersonic except very close to the star. The 
wind material will be spread more equally over all lateral regions. This will avoid the problems which previous 
models |32| had with the unobserved polarization of the stellar light [|6j . 

Previous wind models had to use rather slow rotation and rather strong magnetic fields to produce the high 
observed terminal velocities and to avoid conflicts with the not observed polarization of the stellar light. But 
such models lead to a rapid spin-down. We can now avoid this using a wind model with rapid rotation and 
moderate magnetic fields. In this case we have strong mass loss combined with a sufficiently small angular 
momentum loss. In such a model the mass loss will be concentrated in the equatorial regions. But our bent 
fluxsheets will distribute the wind material in the polar regions as well, so that we can avoid the polarization 
problem. Additionally our outward running waves will reproduce the high observed terminal velocities without 
requiring a stronger magnetic field. These two models in combination help a lot to fit theoretical models for 
magnetic winds to the observations. 

Nevertheless we are not a the end of the development. We are only at the beginning. And there is still a lot 
of interesting work to do. The next, direct step should be a deeper mathematical and numerical analysis of our 
fluxsheet model in order to complete this work. In general both of our models can be improved. For our analysis 
of waves in magnetic winds the next step is to look at the nonlinear case and to demonstrate that our predictions 
for the behavior of magnetic shock fronts are correct. This work could produce quantitative predictions for 
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the X-ray emission, for the nonthermal radio emission, and for the cosmic ray emission. This would improve 
our understanding of the microphysics of these winds. Our fluxsheet model can be improved beyond the steps 
already mentioned here and in Chap. |t] in two directions: (1) A more detailed description of the radiation driving 
mechanism can give a better understanding of optical thick winds [E3L And the non-radial component of the 



radiation force can bend the fluxsheets even more towards the poles 57 . It should therefore be included in our 



model. This might be important in order to show, that the wind material is really equally distributed around 
the star, and we therefore see no polarization. (2) The second major direction in the future development of the 
fluxsheet model is the implementation of a real two dimensional treatment of the MHD equations to get a better 
understanding of the shape of the fluxsheets and of the large scale structure of the wind. 

Another large area of future work is the application of the fluxsheet model to other winds. On purpose we 
did not integrate a specific geometrical shape of the fluxsheets into the fluxsheet model. With few modifications 
it is possible to use the model with a different fluxsheet geometry. This allows e.g. to model the wind from 
accretion disks in young stellar objects or in AGNs. This opens a wide range of future projects. 

Thus at the end of this thesis we have the exciting situation of science: Answer one question, and you find 
three new ones! 



Appendix A 

The equatorial wind as limiting case 



In this appendix we derive the equations for the equatorial wind from Chap. Q as a limiting case of the wind 
equations from Chap. 0. This is an important consistency check for the new model. Additionally it allows to 
relate results calculated with the equations of Chap. || to results calculated with the equations of Chap. 0. The 
equations of Chap. [| can be obtained by using the fluxsheet shape functions 

x(u,rj) = - (A.l) 
u 

6(u,rj) = arccos[a(M)?7], (A. 2) 

where a(u) is the area function for the equatorial fluxsheet we have introduced in Chap. |J. The equatorial 
fluxsheet is specified by 

j] = 0. (A.3) 

The symmetry between northern and southern hemisphere requires that our self-similarity functions Y, F, uj, 
PAc, and BpAc vanish in the equatorial plane. For Y we choose 1 in the equatorial plane. We find then for the 
geometry expressions 
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All other relevant geometrical quantities are zero. We find then for the auxiliary quantities (Eqs. 7.55 , 7.60 - 7.70 ) 
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reproduces our results from chapter |5j 
(Eqs. |5.l| - |5.4 ). It might now appear that the equatorial version of the wind equations is superfluous since we have 
the general version. But they were an important step in the development of the generalized model. Additionally 
they are much simpler, so that there numerical implementation is faster and more robust against numerical 
roundoff errors. Therefore it is useful to start numerical computations with the equatorial wind equations. 



where we have used A = w 2 /""o- Tnis result ( E q s - 7 - 56 , A.26, A.32| & A.37 



Appendix B 

Numerical considerations 



Before we can integrate Eq. 7.56 to find a unique wind solution, we have to specify some numerical parameters: 
(1) We have to describe the geometrical shape of the fluxsheet by specifying r*, x(u, 77), and 6(u, -q). This defines 
|V7j| as well. (2) The next step is to specify the self-similarity functions Y{rj) and pa c (v)- These functions 
have to be guessed for the first calculations. Then we can derive them from previous calculations until we 
get a self-consistent combination of self-similarity functions and wind solutions. (3) The star itself is described 
by its luminosity L*, its mass M, and its equatorial radius i? cq , its (equatorial) magnetic field -B p oEq, and its 
equatorial rotation rate a 10 t,cq. (4) The physics in the wind depends on the CAK parameters a ca k, fccak and on 
the temperature stratification T(u). All models presented in this thesis assume an isothermal wind. (5) Finally 
some of theses quantities may depend on latitude. For all numerical models presented here we ignore differential 
rotation a I0 t(v) = a rot,cq and a lateral dependence of the magnetic field. But the lateral dependence of the 
stellar radius can not be ignored. Some wind models in this thesis have a rather high rotation rate, which will 
lead to an oblate star. Without going deep into the theory of stellar structure we can describe the stellar surface 
as an equipotential surface of 

Q , GM fi 2 r 2 sin 2 

$r,0 = . B.l 

r 2 

The radius R(rf) for the base of the wind is given by the point where the fluxsheets intersects the stellar surface. 
We can find this point by solving 

<Z>(r*x(uo(ri), rj),6{u {rj), r,)) = $ (i? oq , |) (B.2) 
for uq (77). Finally we have to specify eigenvalues as in the equat orial case. We chose again the position of the 



Alfvenic point uac and the approximate mass loss rate (Eq. 7.37 ). Using 
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the wind velocity j/ac and the density pAc at the Alfvenic point can be found 



In the case of a cold wind (v s = 0) we integrate Eq. 7.56 from Alfvenic critical point ua c to the base of the 



wind ito and to infinity. We can obtain ua c by fitting a reasonable initial wind velocity vo at the base of the 
wind. 

In the case of a warm wind (v s > 0) we integrate from the inner critical point to the base of the wind and 
to the Alfvenic critical point and from the outer critical point to the Alfvenic critical point and to infinity. We 
can now fit the initial wind velocity at the base of the wind and a continuous wind acceleration at the Alfvenic 
critical point. These two conditions allow to find uac and the mass loss rate. 

One important point has yet not been taken into account properly. The plot for the transverse component 
of the magnetic force shows a singularity at the Alfvenic point. This is of course not physical. The reason for 



this is the uncompensated (1 — Mp) term in the formulas for the various derivatives of y^ (Eqs. 7.64 



-7.67) 



7.73). The 



These derivatives enter the wind equation (Eq. |7.56| ) and the equation for the transverse forces (Eq, 
wind equation has a singularity at the Alfvenic point. But this is already known from the equatorial theory. 
The coordinates of the Alfvenic point are chosen so that we get a solution with finite y and y$ at the Alfvenic 
point. But nevertheless this is a problem which has to be fixed before really reliable results for winds in bended 
fluxsheets can be calculated. We might derive an additional condition for the wind solution at the Alfvenic point 
from the proper treatment of the transverse Euler equation. 
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Another challenging point is the computation of the inner and outer critical points. The wind equation derived 
in this chapter is much more complicated than the equatorial equations derived in the previous chapters. We 
derive here the equations for the CAK-type critical points in the non-equatorial wind using the same argument 
as in the case of the equatorial wind. At the CAK-type critical points we have 
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With Eq. B.6 y can be eliminated from Eq. B.7 leading to 
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Eqs. B.5, B.6 



B. 



form a set of three nonlinear equations with three unknowns: u, y, and y. We can therefore 
expect at least locally unique solutions. Since these equations are very complicated and essentially non- algebraic, 
it is not possible to find an analytic solution. We use a numerical method based on Newton's algorithm. This 



algorithm requires the Jacobian matrix of the terms on the right hand sides of Eqs. B.5| , B.6 
to it, y, and y. This requires the second derivatives of A, B, and C . 



B.8 with respect 



Since Eq. B.8 is very complicated it seems reasonable to evaluate the derivatives of A and B numerically. 



But we have to consider that the Jacobian matrix, which contains the second derivatives of A, B, and C, is even 
more complicated and requires certainly numerical evaluation. Tests have shown that evaluating the first and 
second derivatives of A, B, and C numerically causes severe numerical roundoff errors, which make it impossible 
to find the critical points. The reason is that numerical derivatives contain the difference of two nearly identical 
quantities. This l eads to a lo ss o f precision. In our case the situation is worse due to the fact that the right 
hand sides of Eqs. B.5, B.6 & B.8 vanish at the critical point as well. I.e. they are differences of nearly identical 



quantities close to the critical points. We control these problems by using quadruple precision for all numerical 
calculations. Additional we use numerical derivatives only for the Jacobian matrix. Nevertheless our code has 
some times problems finding the outer critical point numerically. But it is not clear yet, whether this is an 
numerical or a mathematical problem. See Sect. 7.£ for more details. Finally it should be mentioned that in 
general a forth condition for the inner and outer critical point can be found. A 'perfect' wind solution should 
fulfill the transverse Euler equation everywhere. In this case we can demand that analogous to Eq. B.7 the total 
derivative of the transverse Euler equation with respect to u vanishes. Further details are discussed by Falcke 
|lS| , Sect. 5.3]. Our solutions do not obey the transverse Euler equation in general. Therefore we can not use 
this condition for our critical points. 

The analytical expressions for the first derivatives of A, B, and C are given now. The formulas are rather 
complicated and offer certainly some space for simplifications. 
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